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PREFACE 


The Global Modeling and Simulation Branch (GMSB) of the Laboratory for 
Atmospheres (GLA) is engaged in general circulation modeling studies related 
to global atmospheric and oceanographic research. The research activities 
are organized into two disciplines: Global Weather/Observing Systems and 

Climate/Ocean-Air Interactions. 

During the past year, a very large number of research projects have been 
completed as indicated by the list of refereed papers included at the end of 
this publication. The purpose of the review is to present a brief overview of 
this research, and to document the status of other projects currently underway. 
In addition, abstracts of the 1984 Summer Lecture Series jointly sponsored by 
the Global Modeling and Simulation Branch and the University of Maryland are 
also included. 
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I. GLOBAL WEATHER/OBSERVING SYSTEMS 


A. ANALYSIS AND FORECAST STUDIES 



Analysis of Ocean Surface Wind Fields Using Seasat-A Scatterometer Data 

E. Kalnay and R. Atlas 


In addition to the Seasat forecast impact studies previously conducted 
with the GLA model, a study was also made of the global analyses of surface 
wind stress, wind stress curl, and surface sensible and latent heat fluxes, 
which are generated as a by-product of the assimilation process. Both instan- 
taneous analyses of these quantities and 15-day averages of each of these 
fields were found to be in good qualitative agreement with published fields. 

Examples of Synoptic (Instantaneous) Analyses 

Fig. 1 presents the GLA surface wind streamline analysis for 00 GMT 15 

All the meteorological features described by Peteherych et al. 
(1984), are also present in the GLA analysis. Among these note the two typhoons 
off the coast of Japan, southwesterly monsoon flow on the coast of equatorial 
West Africa and the Indian subcontinent, the cyclones and anticyclones in the 
Northern and Southern Hemispheres, the strong westerlies that drive the Antartic 
circumpolar Ocean current, and the trade winds from the Northern and Southern 
Hemispheres meeting in the ITCZ at about 10°N, and others. In addition, a weak 
cyclone at 45 N, was detected in the GLA analysis and not in the hand analysis 
because it is in a region where no SASS data was available due to the SEASAT 
orbital configuration. 

Fig. 2 presents the instantaneous wind stress field for the same date. The 
most notable feature is the intense stress associated with the enormous Southern 
Hemisphere-cyclone pointed out by Peteherych et al . and centered at 32°S and 135°W. 

15 Day Averaged Surface Wind and Wind Stress Fields 

Fig. 3a, b presents a comparison of the long term mean vector surface 
winds in the GATE area during September adapted by Burpee and Reed (1981), and 
our GLA analysis surface wind for the 15 day period 6-20 September 1978. There 
is very good agreement between these fields, with a col at about 10°N and 27°W. 

West of this col the trade winds converge into the ITCZ which is located at 
about 10 N. East of the col, the winds are towards the east, converging towards 
a broad cyclonic circulation located over the maximum surface temperature at 
about 20°N and 10°W. The cyclone center in our analysis is displaced towards 
bhe Greenwich meridian. The mean (vector) wind speeds are also in reasonable 

with the larger values in the GLA analysis probably due to the short 
averaging period. 

REFERENCES 


Kalnay, E. and R. Atlas, 1985: Global analysis of ocean surface wind and wind 

stress using the GLAS GCM and Seasat scatterometer winds. Accepted by J. 
Geophys. Res. ' 

Burpee, R. W. and R. J. Reed, 1981: Chapter 4, "Synoptic Scale Motions". NSF 

GATE Report. 

Peteherych, S., M. G. Wurtele, P. M. Woiceshyn, D. H. Boggs, and R. Atlas, 1984: 

First global analysis of Seasat scatterometer winds and potential for meteor- 
ological research. Frontiers of Remote Sensing of the Oceans and Troposphere 
from Air and Space Platforms . 575-586. 
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Fig. 3 Mean surface winds in the GATE area during September 
(a) adapted from Burpee and Reed (1981), (b) 15 day 
average from GLA analyses for 6-20 September 1978. 
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Forecast Impact Simulation Studies 


R. Atlas, E. Kalnay, W. E. Baker, J. Susskind, D. Reuter, and M. Halem 
Laboratory for Atmospheres, NASA/Goddard Space Flight Center 

Greenbelt, MD 20771 


A series of simulation experiments is being conducted as a cooperative 
effort between the European Centre for Medium Range' Weather Forecasts (ECMWF) , 
the National Meteorological Center (NMC) and the Goddard Laboratory for Atmos- 
pheres (GLA), to provide a quantitative assessment of the potential impact 
of proposed observing systems on large scale numerical weather prediction. 

For these studies an advanced analysis/forecast simulation system has been 
developed which provides for a more realistic assessment of the impact of 
proposed observing systems than was possible in earlier studies. This system 
consists of four elements: .1) An atmospheric model integration to provide a 

complete record of the "true'' state of the atmosphere (called nature). This 
record is then used to fabricate observational reports and to evaluate ana- 
lyses and forecasts. 2) A conventional data assimilation cycle that is used 
as the "control experiment." The control experiment is like an operational 
analysis/forecast cycle based on conventional observations, except that it 
makes use of fabricated conventional data obtained from the nature run to 
produce the analyzed fields. 3) A satellite data assimilation that differs 
from the control in also including fabricated satellite data incorporated in 
the forecast-analysis cycle. 4) Forecasts produced from both control and 
satellite initial conditions. Comparison of these forecasts with nature 
provides an assessment of the impact of satellite data. 

Previous simulation studies have been characterized by the use of the 
same model to simulate "nature" and observations and to produce forecasts. 

This "identical twin" problem can distort the conclusions derived from such 
studies. In the present study, we attempt to avoid these limitations by de- 
signing a more realistic simulation system and calibrating its results by 
comparison with real data experiments performed with a similar system, and 
by simulating the expected accuracy and characteristics of observational 
systems. In order to avoid the "identical twin" character of previous stud- 
ies, the high resolution (1.875° x 1.875° x 15 levels) ECMWF model is used 
as nature, and the 4° x 5° x 9 level GLAS model is used as the assimilation 
and forecast model. 

For all of the experiments which have been completed to date, the na- 
ture run is a twenty-day integration from 0000 GMT 10 November 1979 using 
the ECMWF model. All types of FGGE and conventional data were simulated by 
NMC by interpolating the nature fields to observation locations and adding 
assumed random or systematic errors to the interpolated values. Only satel- 
lite temperature soundings were assumed to have systematic errors. LIDAR 
wind profiles were simulated at the TIROS observational locations with l-3m 
sec"l accuracy. Wind profiles were not generated at levels below which the 
integrated cloud amount exceeded 90%. 
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Experiments have been conducted to calibrate the simulation system and 
determine its realism, and to begin to assess the relative impact of tempera- 
ture and wind profile observing systems. To this end, two real data assimi- 
lation cycles, a control (which included only conventional data) and FGGE 
(which included conventional data plus satellite temperature soundings, 
cloud-track winds and drifting buoys and were performed for the period 0000 
GMT 10 November to 0000 GMT 25 November 1979. Six simulated data assimila- 
tion cycles, control, FGGE, control plus TIROS temperature soundings, control 
plus perfect temperatures, control plus cloud-track winds, and control plus 
wind profiles were then generated for the same period. 

The NMC analysis for 0000 GMT 10 November was used as initial conditions 
for the real data assimilation cycles. Initial conditions for the simulated 
data assimilations were provided by a real data control assimilation from 
0000 GMT 4 November to 0000 GMT 10 November. Eight five-day forecasts were 
generated from each assimilation at 48 h intervals beginning on 11 November. 

The results from these experiments indicate that 3-dimensional wind 
profile data are more effective than temperature data in controlling analysis 
errors. If the proposed accuracies and coverage for LIDAR wind profiles can 
be achieved, the experiments show that a very significant improvement in 
Southern Hemisphere forecast skill will result. In the Northern Hemisphere, 
the impact is smaller but still significant: a case study showed a major 

improvement in the prediciton of a storm over the United States, which was 
poorly forecasted with the simulated current system. The use of either 
perfect temperature soundings or wind profile data resulted in an improved 
cyclone forecast. 

Further experiments with this system will be used to assess the degree 
of forecast improvement that might result from different configurations of 
Lidar observing systems. Improved methods for assimilating future observa- 
tions will also be studied with the simulation system. 
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The Impact of Scatterometer Data on Limited-Area Model 
Predictions of the QE II Storm 

D. Duffy and R. Atlas 

Experiments have been conducted to assess the impact of Seasat-A scatter- 
ometer (SASS) data on high resolution model predictions of the QEII storm 
and to explore different approaches for utilizing SASS data. For all of 
the experiments, the model used is the limited-area model described by Duffy 
(1981) with a horizontal resolution of 100 km. 

Numerical weather forecasts were made from initial conditions which included 
and excluded SASS data. The initial condition for the control run was the 
global analysis from the National Meteorological Center (NMC) which we had archived 
at the resolution of 4° latitude and 5° longitude. This field was bilinearly 
interpolated to the model's 100 km by 100 km grid. For the SASS experiments, 
the 1000 mb winds were replaced by subjectively dealiased SASS winds wherever 
they were available. Furthermore, SASS winds were allowed to influence the 
upper-level winds through the use of a vertical correlation function (Yu and 
McPherson (1984)). If the upper-level winds were not modified, the effects 
from the upper levels soon overwhelmed the corrections that we had made at the 
surface. 

In Fig. 1, we present the sea-level pressure and 1000-500 mb thickness at 
the initial time of the experiments, 12 GMT 9 September 1978. We were primarily 
interested in the evolution of the nascent low located just off Nova Scotia. 

The results from the two different numerical forecasts are presented in 
Figs. 2 and 3. In both cases the low located just off Nova Scotia has moved 
out into the open waters of the North Atlantic Ocean. In the case of the 
control, the low has intensified only slightly. The central pressure contour 
at 36 hrs is 1000 mb with peak winds of 21 m/sec. These results are slightly 
better than those given by the operational limited-area, fine-mesh model run 
at the NMC and the Navy's operational model. 

From Fig. 3 we see that the inclusion of SEASAT data has had a dramatic 
effect on the cyclogenesis. The central pressure contour at 36 h is now 988 mb 
and the peak winds are 37 m/sec. A study of the results indicates that the 
deepening of the low occurs because of larger amounts of convective precipita- 
tion in the SEASAT experiment than the control. 

In Fig. 4 we show the observed NMC hand analysis of the QEII, valid at 00 
GMT 11 September 1978. From this figure we see that the forecast with SEASAT 
data is closer to the observed storm than the control. However, there is still 
a considerable position error. 

These results show that scatterometer data can be important in defining 
the initial conditions that result in rapid maritime cyclogenesis. Currently, 
wbrk is under way to determine improved ways to couple the surface winds to the 
upper levels in global atmospheric models. 
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Objective Determination of Heat Wave Patterns 


N. Wolfson and R. Atlas 


A simple index was developed in order to objectively characterize a heat 
wave pattern in either analyses or forecasts of 500 mb height. 500 mb height 
was chosen because this parameter is usually among the most accurate fields 
predicted by numerical models. In addition, it is highly correlated with 
the surface air temperature on many time scales. 

The heat index IN is designed to be applied to any region of interest and 
to be used for any objectively analyzed data and specifically those of NWP 
models* IN is given by 

IN = 77j7 I / H 500 " H 500 ds dt (1) 

b T S 

where H 500 is the 500 mb height 

h 500 the 500 mb climatological value of H 500 

Z, t are location and time 

ds is an area element (represented by a grid 
box for an NWP model) 

S is the total area of interest 
T is the period under study 

The index can be calculated_f or any time period, day, month or season 
with usage of the corresponding HsqoU). It can be used to study results 
of NWP models by substituting forecasted values of H 50 o(^,t) into (1) or to 
diagnose observations by using objectively analyzed values of H 5 oo(£,t). 

In formulating the index, it was assumed that point values of height are 
representative of field values in an area element. The height field spatial 
autocorrelations at 500 mb are very high for distances of - 250 km (Buell, 
1972). Thus for the distances involved here, within an area element, point 
values can be considered to be highly representative for the areas involved. 

The first integration over the area transfers the local anomalies to an 
areal measure. It enables the determination of whether a distinct geographical 
region is experiencing a heat wave or whether there is simply a local event. 

The second integration in time assesses the temporal average anomaly of the 
event . 

To evaluate the seasonal represent iveness of the index for the summer 
season, June-September seasonal values for the years 1963-1977 were correlated 
with average seasonal temperatures for the northern and southern Plains of the 
U.S. The index was calculated from the U.S. National Meteorological Center 
(NMC) 500 mb analysis after interpolation to the GLA grid. Average seasonal 
temperatures were calculated by averaging monthly average temperature which 
were taken from U.S. Dept, of Commerce (1981). The correlation coefficients 
were .85 and .82 for the North and South Plains respectively. 
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Fi.§ure 1, which presents the daily variation of the heat wave index for 
the North and South Plains, clearly shows the different nature of the heat wave 
in these two regions. The horizontal lines are the seasonal threshold values, 
i.e. higher values indicate a heat wave. The heat wave in the southern plains 
was rather continuous with short breaks, the northern plains experienced numerous 
spells of hot weather which were shorter than the southern spells. 

As the area integral in (1) is reduced to one grid point one obtains the 
measure for the local heat wave. A grid point by grid point calculation of 
monthly values enables the monitoring of the local monthly variability. In 
order to assess the severity of the heat wave at a local grid point £ a "t” 
measure has been applied thus 

tOO = - IN (*> (2) 

®eU) 

where IN(£) is the monthly index at location 9 . 

og is the monthly standard error of IN. 

IN and eg are calculated from the same data, which can be observations 
or forecasts. As both are in meters t(£) is dimensionless. 

The numerator in (2) gives the temporal average anomaly. This information 
is not enough to decide whether the heat wave was characterized by a series of 
events with high anomalous values with frequent in between breaks or by a 
persistant anomaly with no breaks. Dividing by the standard error allows us 
to give weight to the heat wave character (small standard errors are indicative 
of persistant H(£,t) - H(£,t) values during the time perids under consideration). 

This concept is close to the M t” student test checking the null hypothesis 
that IN( £) equals zero. Interpretations of t(£) as "t” student test values is 
very complicated due to the difficulties in determining the effective number of 
degrees of freedom and other statistical constraints. However, it is obvious 
that higher t(£) values are indicative of a stronger heat wave as a combination 
of the two factors of intensity and persistence. 




Fig. 1. Daily changes of the heat wave index in the North and South Plains. 
Threshold values for the heat wave are shown as horizontal lines 
(gaps indicate missing data). 
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Wave Structure Associated with the Summer 1980 Heat Wave and Drought 


N. Wolfson and R. Atlas 


A study has been made of the wave structure associated with the summer 1980 
heat wave and drought in relation to other years. This study was performed for 
the 500 mb level for the months of June-Sept. 1963-1977 and 1980. The data 
was obtained from the National Meteorological Center (NMC) final analysis inter- 
polated to the GLA 4° latitude x 5° longitude grid. For each year the 500 mb 
height was spectrally decomposed. Then waves 1-3 and 4-7 were combined. The 
decompositions were carried out for latitudes 30°N to 54°N every 4°. Fig. la 
presents the longitudinal location of the ridge lines for waves 1-3 for the 
period June-Sept. of 1963-1977, 1967 and 1980. One may see that there are three 
ridge lines, one over Tibet and central Asia and the other two over the central 
U.S. and eastern Atlantic. The ridge over the U.S. is very stable and does 
not change much from year to year. This is indicated by the bar line on Fig. la 
which indicates one standard deviation in the yearly location. The deviation 
from the mean is usually smaller than 20° longitude. 

As indicated by Fig. la the summer of 1980 has a different pattern and 
shows large deviations from climatology. The ridge over the central U.S. 
north of 38° shifted by 10°-20° eastward. This shift is more apparent in the 
northern latitudes and is larger than the standard deviation. 

The ridge lines of waves 4-7 are presented in Fig. lb. Ridges are observed 
near main mountain chains i.e. the Alps, Caucasus, Tibet and the Rockies, 
indicating the importance of these ranges in influencing the summer circulation. 
One may note that over the central U.S. the planetary waves 1-3 and the shorter 
4-7 waves support each other, with both contributing to the standing wave. 

The ridge line follows the eastward slope of the longitude of the Continental 
Divide supporting Wallace's (1983) hypothesis of the effect of orography on 
the summer circulation. It is evident that summer 1980 was highly anomalous 
in that the ridge north of 35°N shifted eastward by ^15° longitude, while the 
climatological standard deviation is only 5°-12° longitude. 

Further evidence of the relationship between the location of the planetary 
waves and average seasonal temperature is gained by analyzing the wave structure 
for 1967, which was the coldest summer of the 15 years analyzed in the northern 
plains. North of 38N, waves 1-3 were shifted westward (Fig. la) by 15°-25° 
longitude, while waves 4-7 showed a 10° average shift from the climatological 
location in the same direction. Furthermore, the 1967 summer was the warmest 
year on record between 1963-1977 in the pacific region. 

In order to study the relative importance of the 1-3 and 4-7 waves during 
the heat wave episode we have plotted their amplitude for the period June 
through September 1980. 

95°W is a representative longitude for the Great Plains. At 30°N,95°W the 
main contribution of the heat wave is from waves 1-3 whose amplitude is almost 
always larger than those of waves 4-7. However, both of these wave bands had 
positive amplitudes and supported each other during this period. This is similar 

at 34 °N but the contribution of waves 4-7 was more important for the maintenance 
during June and the second part of August. At 38 °N tne combination of 1-7 is 
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Preliminary Results from Numerical Experiments on the 
Summer 1980 Heat Wave and Drought 


N. Wolfson, R. Atlas, and Y. Sud 


During the summer of 1980, a prolonged heat wave and drought affected the 
United States. The heat wave began in the southern Great Plains in June extendec 
northward and westward in July and then in August and September weakened and 
retreated to the Southeast. 

A preliminary set of experiments has been conducted to study the effect 
of varying boundary conditions on the GLA model simulation of the heat wave. 

Five 10-day numerical integrations with three different specifications of 
boundary conditions were carried out: a ,, control ,, experiment which utilized 

climatological boundary conditions, an "SST" experiment which utilized summ- 
er 1980 sea-surface temperatures in the North Pacific but climatological values 
elsewhere, and a "Soil Moisture" experiment which utilized the values of Mintz- 
Serafini for the summer 1980. The starting dates for the five forecasts were 
11 June, 7 July, 21 July, 22 August, and 6 September of 1980. These dates 
were specifically chosen as days when a heat wave was already established in 
order to investigate the effect of soil moistures or North Pacific sea-surface 
temperatures on the model’s ability to maintain the heat wave pattern. The 
experiments were evaluated in terms of the heat wave index (Wolfson and Atlas, 
1985) for the South Plains, North Plains, Great Plains and the entire U.S. In 
addition a subjective comparison of map patterns has been performed. 

Table 1 presents the average ten day index values for the four regions. 

The differences between the various experiments are very small at the beginning 
of the integration (not shown), and usually develop after the fourth or fifth 
day. One may note that the SST influence has been usually to decrease the index 
value. This might be due to the cold ocean temperatures prevailing in 1980. 

One may place confidence in these SST fields only for short-range forecasting 
because a longer period of integration would demand global SST distribution to 
be introduced. On the other hand, introduction of 1980 real soil moisture values 
over the U.S. has enhanced the heat wave. This effect is much more prominent 
over the northern plains than over the southern plains. 

Table 2 summarized the impact of boundary conditions on the heat wave index. 
In addition, "t" student test values are shown to indicate the significance of 
the results. 

Additional model experiments are being conducted to study the initiation 
and breakdown of the heat wave and the effect of global sea-surface temperatures. 

REFERENCE 

Wolfson, N. and R. Atlas, 1985: Objective determination of heat wave patterns. 

(This issue). 
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Table 1: Average ten day index values 


Experiment 

Control 

SST 

Soil Moisture 


South 

Plains 


6/11 

56.06 

52.38 

55.16 

7/07 

58.8 

57.69 

62.38 

7/21 

46.7 

45.55 

44.56 

8/22 

64.52 

55.72 

58.88 

9/06 

60.70 

58.33 

64.23 


North 

Plains 


6/04 

3.59 

-1.68 

4.50 

7/07 

50.8 

47.18 

54.30 

7/21 

50.99 

52.26 

54.15 

8/22 

47.0 

44.18 

53.71 

9/06 

54.0 

47.76 

58.81 


Great 

Plains 


6/11 

26.2 

21.69 

27.54 

7/07 

54.3 

51.73 

57.79 

7/21 

49.1 

49.36 

50.00 

8/22 

54.58 

49.17 

55.95 

9/06 

56.90 

52.33 

58.93 


U. 

,s. 


6/11 

9.20 

7.56 

9744 

7/07 

10.60 

13/21 

14.08 

7/21 

43.01 

43.49 

42.08 

8/22 

17.71 

14.80 

20.79 

9/06 

40.08 

38.72 

44.63 


Table 2: Impact of boundary conditions on heat wave index 

(significant values are underlined) with “t” value in lower 
right side. 


Region 


South 

Plains 

North 

Plains 

Great 

Plains 


Average 
control value 

57.36 


41.28 


48.22 


Impact 
of SST 


-3.42 


2.68 


-3.34 


3.30 


-3.36 


3.71 


Impact of 
Soil Moisture 


.31 


1.17 


3.42 


4.05 


1.83 


4.54 

































B. SATELLITE OBSERVING SYSTEMS 




Improved Cloud and Surface Fields Derived from HIRS2/MSU 
J. Susskind and D. Reuter 


The methods to retrieve cloud parameters from HIRS2/MSU data, described in 
Susskind et al. , (1984), have been shown to produce qualitatively reasonable cloud 
fields but to underestimate cloud amount. Improvements have been made to the 
cloud algorithm which allow for better determination of low clouds. These 
improvements involve retrieving cloud fields at a higher spatial resolution, 
which decreases the problems arising from multiple layer clouds, and increasing 
the number of channels used to determine the cloud parameters from 2 15 ym lower 
tropospheric sounding channels to 5 channels, including the 11 ym window and 2 
15 ym mid and upper tropospheric sounding channels. The use of the 11 ym window 
channel was very important in identifying low level clouds. 

Another way of checking for the existence of low clouds, at least over the 
open ocean during the day, is by looking at the reflected solar radiation in HIRS2 
channel 20, sensing the near IR and red regions. Channel 20 was not used in the 
cloud algorithm because of the desirability to use the same algorithms day and 
night in order to generate meaningful statistics of day-night cloud differences. 

Global one day cloud fractions and cloud top pressure fields have been 
derived from the HIRS2/MSU data. These fields are completely consistent with 
visible and infra-red imagery taken the same time from AVHRR but contain consider- 
ably more information than the AVHRR images. In addition to cloud height informa- 
tion, which is not apparent from the imagery, there is also good cloud information 
at high latitudes where the visible imagery on AVHRR sees either no signal (per- 
petual night) or shows reduced contrast because of background ice and snow cover, 
and the II ym imagery on AVHRR is hampered because of the low thermal contrast 
between the clouds, which are usually low, and the cold surface. 

Monthly mean fields of cloud fraction, cloud top pressure, and cloud top 
temperature, for both day clouds (3 PM) and night clouds (3 AM), and their 
differences, have been prepared for a number of months. Average monthly mean 
cloud fractions obtained in this fashion are about 40%, with 3 AM fields being 
about 5% more cloudy than 3 PM fields. These cloud fractions are effective, 
because the clouds have have been treated as black. They represent the product 
of the cloud fraction and the cloud emissivity, which averages about .7 in the 
15 ym region (Paltridge and Platt, 1976). 

Monthly mean "albedo” (actually per cent reflected solar radiation at a 
specific angle in a specific wavelength interval) fields have also been prepared 
by averaging the observations in HIRS channel 20, divided by the cosine of the 
solar zenith angle. These fields show excellent agreement with the day-time 
cloud fields, indicating most clouds which are apparent in visible imagery are 
also found by the infra-red algorithms. The albedo fields also show features 
attributable to deserts and ice and snow, though it is difficult to distinguish 
the ice and snow from the cloudy regions around them from the visible imagery 
directly. 

The cloud algorithm, which does not use the visible data, provides the ability 
to select scenes which are clear or cloudy. Monthly mean fields of "cloud free 
albedo" were generated by eliminating all scenes from the instanteous albedo field 
considered to be cloudy and averaging the albedo only for cloud free scenes. 
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This procedure serves two purposes. First, it checks the extent that the cloud 
algorithm may have missed clouds apparent in the visible, and second, if success- 
ful, it provides a reasonable estimate of surface albedo once atmospheric correc- 
tions can be made. The "cloud free albedo" fields indeed showed only slight trace 
of low cloud contamination, whose sources we have identified and are correcting. 

In addition, desert features were readily apparent and clear indications of ice 
and snow fields emerged from the cloudy mid-latitude and polar regions. These 
ice and snow fields were quite consistent with those derived from the surface- 
emissivity at 50.3 GHz (channel 1 of the MSU, which has a surface transmittance 
of .7) and the retrieved ground temperature. The agreement of these ice and snow 
fields, neither of which could have been derived by single spectral region sound- 
ing, shows the power of a combined infra-red, microwave, and visible monitoring 
system. 

Another extremely important parameter derived from the HIRS2/MSU soundings 
is the monthly mean difference of the 3 PM and 3 AM ground temperatures. This 
difference is related to the forcing function, given by incoming solar radiation, 
and the response function, given by the thermal inertia of the ground, which 
depends on soil moisture and evapotranspiration rate. Clouds affect the day- 
night temperature difference in two ways. The direct effect of clouds is to 
block incoming solar radiation during day, which slows daytime warming, and 
reflect upwelling longwave radiation back to the surface at night, slowing 
nocturnal cooling. The indirect effect of clouds is to cause rain, which 
increases soil mosture and again decreases the diurnal surface temperature 
variation. Indeed, the day-night surface temperature difference fields are 
quite consistent with the retrieved cloud cover fields for all the months 
processed, especially in the regions of convective cloudiness. 

The effects of clouds on incoming and outgoing radiation are, to first 
order, not a factor in the day-night ground temperature fields because retrievals 
are not performed under full overcast conditions and retrievals are performed 
only in the clear portions of partially cloudy scenes as a result of the cloud 
correction algorithm. The day-night temperature difference fields then should 
be a measure of rainfall and soil moisture. Mintz et al . (1985) have developed 
an energy balance model relating the monthly mean day— night ground temperature 
field derived from HIRS2/MSU to soil miosture and evapotranspiration and have 
derived reasonable monthly mean fields for these parameters. 

In summary, HIRS2/MSU data provides the possiblity of obtaining accurate 
fields of a number of important climate parameters, all derived together in a 
consistent fashion. This system, launched at the end of 1978, will be in place 
through the 80 f s. The opportunity exists to develop a 10 year climate data set 
from similar instrumentation. Our immediate plans are to process all of 1979 
with a slightly improved version of the current system and possibly process all 
data up to 1984. 
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Intercoraparison of Sea-surface Temperature Anomaly Fields 
J. Susskind and D. Reuter 


Sea-surface temperature anomaly fields produced from analysis of HIRS2/MSU 
data at GLAS for the months November 1979, December 1981, March 1982, and 
July 1982 were compared to those produced from data from ships, AVHRR, and SMMR 

in a series of NASA sponsored sea-surface temperature workshops. The methods 

used by GLAS to analyze the HIRS2/MSU data were similar to those in Susskind 
et al. , (1984), but improvements were made to the algorithm used to retrieve sea 
and land surface temperatures and in the cloud filtering algorithm. One of the 
improvements involved generating weights for a given sounding, which decreased 
with increasing cloudiness and increasing humidity. These weights should be 
used in generating monthly mean fields. The workshop generated statistics 
for the weighted HIRS2/MSU field only for the month of July 1982. 

Statistics generated by the workshop for the last two months, in which the 

latest improved algorithm (called HIRS version 2 in the workshop) was used to 
process the HIRS/MSU data, are shown in Tables 1 and 2. Climatology is also 
treated as a sensor for comparison. The values given are C, the correlation 
coefficient, B, the bias (sensor-ship), S, standard deviation, and N, the number 
of grid points sampled (2° x 2° grid) comparing fields from a number of sensors 
with ship fields. Statistics are given for global colocations and also for 
those in the North Pacific and North Atlantic Oceans. Results are given for the 
unsmoothed field, and also for 9 point smoothed fields, in which case, both the 
ship field and retrieved field were smoothed and compared to each other. 

To interpret the significance of the standard deviation from ships, an 
anomaly field produced by a sensor can be considered skillful if its standard 
deviation from ships is at least as low as that of climatology, which is represen- 
tative of the ocean signal. A better indicator of skill is the correlation of 
the sensor anomaly field with that of ships. Climatology, by definition, has 
zero correlation and zero skill. 

As seen in Tables 1 and 2, smoothing has reduced the noise and increased the 
correlation coefficients for all cases. Part of this is due to random noise 
reduction in the ship field and part due to noise reduction in the field being 
compared to ships. The AVHRR field in March 1982 shows marginal skill in the 
standard deviation sense and good skill in the correlation sense in the unsraoothed 
fields. Statistics for the unsmoothed March HIRS field were not generated at the 
workshop. 

In the smoothed fields, AVHRR shows increased skill in both categories. The 
statistics for the smoothed HIRS fields are almost as good as those of AVHRR with 
regard to standard deviation but the correlations are considerably lower. It 
should be noted that these statistics are for the unweighted HIRS 2 field. The 
weighted field provides better visual agreement with the ship field but statistics 
were not computed in the workshop using the weighted field. The SMMR field is 
much more noisy but shows some correlation skill in the North Pacific Ocean. An 
additional field produced using both SMMR data and ship data shows good agreement 
in the vicinity of ships as expected. Surprisingly, the statistics are only 
marginally better than those of AVHRR and HIRS which did not have the benefit of 
including the ship data used to verify the solutions. 
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When one looks at the actual anomaly fields for March 1982 (not shown), it 
is apparent that March 1982 was indeed very close to climatology, expecially in 
the North Atlantic Ocean. In the North Pacific Ocean, there are small cold 
anomalies centered at about 160° W, 35° N and 175° E, 28° N and a small warm 
anomaly about 115° W, 20° N. The weighted HIRS2/MSU field likewise shows small 
anomalies in the North Pacific ocean, of the correct sign, centered at the 
appropriate locations. The small anomalies in the Atlantic Ocean are in almost 
perfect agreement with those in the ship field. The AVHRR field, on the other 
hand, is shown to indicate extensive areas of cold anomaly in both the North 
Atlantic and Pacific Ocean, sometimes exceeding 1.5° in magnitude. This feature 
is consistent with the cold bias of .3° to .5°C indicated in Table 1 for AVHRR. 

The AVHRR field also has a very large spurious warm anomaly in the region along 
the equator from 75° E to 155° E. Agreement of all fields in the southern 
hemisphere is reasonable, but the AVHRR warm anomalies are larger than those 
indicated by HIRS2/MSU or the ships. The unweighted HIRS field is somewhat warmer 
than the weighted field and has essentially no anomalies in the North Pacific 
Ocean. It is not surprising that small correlation coefficients were found in 
the North Pacific Ocean for the unweighted HIRS Version 2 field. 

Statistics for July 1982 are shown in Table 2. For this month, statistics 
were generated for both the weighted and unweighted HIRS version 2 fields. As 
observed from looking at the statistics comparing climatology with ships, July 
1982 was much more anomalous than other months in the workship with a standard 
deviation of about .6°C about a cold bias (ships-climatology) of .96°C in the 
North Pacific Ocean and .4°C in the North Atlantic Ocean. Therefore signals 
are large and one expects larger correlations between retrieved anomaly fields and 
ship fields. Correlations for AVHRR are similar to those of March 1982 but the 
standard deviations are considerably larger than in the other months. HIRS 
statistics show improved correlation over March 1982, as expected, with standard 
deviations slightly degraded over those of March 1982. Statistics for the weighted 
HIRS field, relative to the unweighted field, show the soundings have become 
colder on the average by about .35°. Results have improved somewhat in the 
Atlantic Ocean and degraded somewhat in the Pacific Ocean. SMMR results have 
improved over previous months, both in the standard deviation sense and in the 
sense of error correlation. The SMMR/ Ship field shows good agreement with the 
ships used to produce the field and shows more of an improvement with regard to 
standard deviation over the HIRS 2 field, which did not have the benefit of 
including ships, than in March. 

The detailed anomaly fields for this month show that HIRS again gives very 
good agreement with the ship anomaly fields, which this time are quite large and 
extensive in both the North Atlantic and North Pacific Oceans. The HIRS field 
does have a few areas of spurious cold anomaly in July, located mainly within a 
thin latitude band running from about 15° N - 25° N over most of the ocean, and 
off the west coast of South American from 15° S to about the equator. These 
spurious anomalies, which are mostly of the order .5° - 1°C, may well be due to 
the effects of aerosols put into the atmosphere by the El Chichon eruption. The 
errors are amplified in the weighted field, which is in general colder. The 
effects of this eruption on derived sea-surface temperatures are much more evident 
in the AVHRR anomaly pattern which shows a large block of spurious cold anomaly 
from 10° N - 35° N running across the entire oceanic area, with magnitudes of 
the order of 2°C. 

In summary, the HIRS2 sea-surface temperature fields are statistically 
comparable to those produced from AVHRR, but the errors tend to be spatially 
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random while the AVHRR errors are more spatially coherent. The effect of the 
El Chichon eruption on derived HIRS sea-surface temperatures appear to be consid- 
erably less than on those retrieved from AVHRR. 


Table 1 

Comparison of Monthly Mean Anomaly Fields 
With Ships > 5/Cell 
March 1982 



I UNSMOOTHED 1 

| SMOOTHED * 



N. PAC. 

N. ATL. 

CLOBAL 


N. ATL. 

c 

mm 

.52 

.67 

.77 

mm 

.80 

AVHRR B 


-.50 

-.29 

-.44 


-.30 

S 


.48 

.42 

.29 

■ 

.21 

N 

hi 

434 

267 

368 

■ 

153 

C 




.55 

.30 

.66 

HIRS2 B 




.29 

.36 

.21 





.31 

.34 

.22 

HUH 




368 

210 

153 

c 

.24 

.37 

-.05 

.15 

.54 

-.09 

SMMR B 

-.21 

.05 

-.76 

-.17 

.13 

-.77 

NICHT S 

1.11 

.99 

1.19 

.79 

.67 

.69 

N 

690 

392 

213 

300 

200 

95 

C 

.58 

.57 


.75 

.75 


SMMR/SHIP B 

.04 

.03 


.07 

.07 


COMPOSITE S 

.47 

.46 


.25 

.25 


N 

438 

394 


207 

203 


C 

.40 


.50 

.79 


.79 

VAS B 

.90 


.89 

.91 


.91 

S 

.56 


.52 

.26 


.26 

N 

109 

1 

106 

51 


51 

C 

.39 

.38 


.70 

.70 


XBT B 

-.27 

1 -.29 


-.47 

-.47 


S 

.89 

.91 


.35 

.35 


N 

242 

! 227 


18 

18 



.00 


.00 

mm 

.00 

.00 

CLIM B 

.09 


-.05 

mStm 

.29 

-.10 

S 

.52 


.42 


.32 

.27 

1 H 

795 

Bfl 

1 267 

mm 

210 

153 


Table 2 

Comparison of Monthly Mean Anomaly Fields 
With Ships > 5/Cell 
July 1982 



UNSMOOTHED 

| SMOOTHED 


GLOBAL 

N. PAC. 

N. ATL. 

CL0BAL 

N. PAC. 

N. ATL. 

C 

.62 

.59 

.66 

.70 

.63 

.84 

AVHRR B 

-.48 

-.37 

-.57 

-.35 

-.17 

-.48 

S 

.79 

.93 

.60 

.52 

.62 

.37 

N 

644 

320 

258 

274 

117 

157 

C 

.49 

.43 

.51 

.78 

.52 

.85 

HIRS2 B 

-.07 

.01 

-.08 

.09 

.14 

.04 

VERSION 2 S 

.69 

.72 

.62 

.38 

.39 

.36 

N 

662 

337 

259 

327 

170 

157 

HIRS2 C 

.52 

.45 

.56 

.79 

.45 

.88 

VERSION 2 B 

-.37 

-.31 

-.37 

-.25 

-.23 

-.27 

WEIGHTED S 

.69 

.74 

.61 

.38 

• 43 

.30 

N 

662 

337 

259 

327 

170 

157 

■H 

.46 

mm 

.32 

.55 

.62 

.66 

SMMR B 

-.43 

mSSm 

-.88 

-.69 

-.39 

-1.07 

NICHT S 

.97 

■H 

.93 

.60 

.48 

.51 


522 

278 

193 

230 

127 

103 

C 

mm 

.74 


.86 

.86 


SMMR/SHIP B 

IB 

-.06 


-.10 

-.10 


S 

■Hi 

.54 


.27 

.27 


N 

316 

282 


137 

137 



.49 


.46 

.42 


.42 

VAS B 

.48 


.50 

.55 


.55 

S 

.46 


.45 

.22 


.22 


92 


88 

38 


38 


.58 

.55 





XBT B 

-.22 

-.23 





S 

.94 

.96 






154 

146 





c 


.00 

.00 

.00 

.00 

.00 

CLIM B 


.67 

.26 

.70 

.96 

.40 

S 


.73 

.69 

.63 

.50 

.64 

N 


338 

259 

336 

179 

157 
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Experiments on the Assimilation of Patches of Data 
Using a Barotropic Model 

F.H.M. Semazzi 


Experiments are performed to examine the impact of adopting various 
procedures in the assimilation of patches of data during the progress of a 
numerical forecast. A barotropic model described by Takacs and Balgovind, 

1983 is used for the experiments. The model is global and its discretization 
is similar to the GLA Fourth Order GCM. The effect of orography is included 
in the model. 

The experiments are of an "identical twin" kind. The "nature" run starts 
from global fields of wind (u and v) and height (h) at 300 mb. They are real 
data of January 1979 00Z which were interpolated from the adjacent sigma levels 
of the GLA GCM and then balanced by the bounded derivative method (Browning 
et al ., 1980, Semazzi, 1985). From this initial state a 6 day forecast is 
performed. At intervals of 3 hours this information is saved over sectors of 
60° width in longitude, from the south pole to the north pole, to simulate 
data coverage and observational frequency of a polar orbiting satellite. 

The initial wind field of the assimilation experiment is constructed from 
simple zonal averages of the u-component of the initial state for the "nature" 
run. The corresponding zonally averaged height field is determined from the 
geostrophic relationship. Every 3 hours the forecast in the assimilation exper- 
iment is interrupted and data is inserted over the appropriate sector region 
from the "nature" run. 

The evaluation of the assimilation proccess is based on the global RMS 
of the difference in height between the "nature" run and the assimilation 
forecast. Fig. 1 displays the limiting case when all the three prognostic 
variables, namely, u,v and h are simultaneously inserted. The two runs become 
identical in about a day when the earth makes a complete rotation. Fig. 2 
shows that the direct insertion of heights has only a small impact. Fig. 3 
indicates that the insertion of winds only is far more effective than the 
height (Fig. 2). Finally, Fig. 4 displays the impact of making geostrophic 
wind correction in addition to the insertion of height* The improvement over 
the results shown in Fig. 2 is dramatic. 

We are presently exploring the prospect of accelerating the assimilation 
by using a combination of the bounded derivative method and the estimation 
theory of stochastic-dynamic systems suggested in Bengtsson et al. (1981). 
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High-Latitude Filtering in a Global Grid-Point Model 
Using Model Normal Modes 

L. L. Takacs, I. M. Navon, and E. Kalnay 


1. Introduction 


The aim of high-latitude filtering in the vicinity of the poles is to 
avoid the excessively short time steps imposed on an explicit time-differencing 
scheme by linear stability due to fast moving inertia-gravity waves near the 
poles. In this paper we apply the model normal mode expansion toward the 
problem of high-latitude filtering in a global shallow water model using the 
same philosophy as that used by Daley (1980) for the problem for large timesteps 
in P.E. models with explicit time integration schemes. 

2. Computational Procedure 

Starting from a non-linear normal mode initialized state, we transformed 
forward at each timestep n the dependent variables (u, v, h) n given to us by the 
model, thus producing a normal mode coefficient set. The model used a conven- 
tional explicit Matsuno time-integration scheme. To obtain the projection 
onto the fast modes, we simply back transformed only those modes whose eigen- 
frequency was greater than the maximum allowed derived from timestep considera- 
tions. This produced Au, Av, and Ah, which were then subtracted from the current 
fields to produce "filtered" new values of those variables. A new integration 
step was then carried out. 

3 . Experimental results 

Figs. 1 and 2 show difference plots after 3 days of the height and wind field 
between three filtering techniques and a control run using no high-latitude filter 
and a small timestep (60 seconds). The three filtering techniques are 1) the norm 
mode filter (NMF) , 2) the prognostic field filter (UVH), and 3) the time-tendency 
filter (DOT). The control run is considered as "truth” in these experiments. 

From these figures we see that the NMF filter has the best simulation of 
the control run, producing very little differences in high latitudes. The DOT 
filter is not much worse. The UVH filter, however, shows major differences in 
high latitudes and begins to spread its influence to lower latitudes as well. 

4) Conclusions 


A methodology has been developed in which fast-moving inertia-gravity waves 
are transformed to normal mode space in order to perform a "surgical" high— 
filter required in global grid— point models. In contrast to other 
high latitude filtering techniques, the new approach leaves the model's Rossby 
modes intact, and is shown to best simulate a control, no-filter run. 
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On the Improvement of Satellite Temperature Retrievals by 
Means of Boundary Layer Models 

R. Boers 


A simulation study was conducted to determine whether boundary layer models 
could be used to improve the accuracy of the satellite temperature retrievals 
close to the earth's surface. We used the simplest type of boundary layer models 
that only consider time changes in the mean boundary layer quantities such as 
equivalent potential temperature and specific humidity, but neglect the details 
of the turbulent structure in the interior of the layer. They can be run with 
relatively few externally specified parameters: surface wind speed and divergence 
rate. Extensive theoretical development of boundary layer models can be found 
in the literature over the last 20 years (Tennekes, 1973, Stage and Businger, 

1981 and many others). Model verification studies (Boers et al . , 1984 and 
others) show that under certain conditions the mean atmospheric structure can 
be simulated very well. 

The following process was carried out to study the impact of boundary 
layer models on the retrieval accuracy. 

a) A set of randomly distributed wind speed and divergence rates was 
assigned to a randomly distributed set of oceanic temperature and 
humidity profiles. 

b) Boundary layer model runs (Stage et al . , 1981) were performed on these 
profiles until a stationary solution was found in boundary layer height 
and temperature. These modified profiles were called the truth. 

c) Brightness temperature (AMTS, HIRS) were computed for the truth and 
temperature and humidity retrievals were performed using these bright- 
ness temperatures. 

d) Boundary layer model runs were performed on the retrievals until a 
stationary solution of boundary layer parameters was obtained. The 
solution was called the modified retrieval. 

e) RMS errors of the retrievals and the modified retrieval with respect to 
the truth were computed. 

Tests included simulations whereby all externally specified parameters at 
step (d) were left the same as those used in generating the "truth", and simu- 
lations whereby each one was perturbed to study the influence of the uncertainty 
of those parameters on the mixed layer solution. Early sensitivity tests 
indicated that in order to acquire meaningful results an accurate scheme was 
needed to compute the radiative sky temperature in case clouds developed at 
the top of the boundary layer. We used the scheme of Harshvardhan et al . (1984) 
to obtain the radiative sky temperature for both the truth and the modified 
retrievals. We can draw the following conclusion from this study: 

a) If wind speed and divergence rate are known exactly for computation of 
the modified retrievals, the modified retrievals are more accurate than 
the retrievals in representing the truth by .2° C to .8° C in the 
lowest 200 mb. 
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b) Typical perturbations of the wind speed by 30-60% do not significantly 
affect these results. 

c) Typical perturbations of divergence by 50-100% decrease the accuracy 
of the modified retrieval to that of the retrieval. 

d) The combined effect of b and c is of no improvement of the modified 
retrieval over the retrieval. 

e) We can significantly improve on the disappointing results of b and c 
by selecting those cases where the truth and the modified retrieval 
both contain clouds, or both contain no clouds, excluding all other 
results. For these cases the modified retrieval is more accurate 
than the retrieval by .4 to .6° C in the lowest 200 mb. 

To illustrate some of our results we consider Figure 1. The four lines 
represent RMS errors in temperature of 1 HIRS-retrieval , 2 HIRS-modif ied 
retrieval, 3 AMTS-retrieval, 4 AMTS-modified retrieval. Excluded are all non- 
matching cloud - no cloud pairs in the truth and the modified retrieval. 
RMS-errors of retrievals and modified retrievals are between .8 and 1.8°C. 

AMTS retrievals are more accurate than the HIRS-retrievals by .2°C except close 
to the surface. The modified retrievals for both HIRS and AMTS improve on the 
retrievals by up to ,8°C close to the surface. 
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Applications of Fuzzy Clustering Techniques to Stratified by Tropopause 

MSU Temperature Retrievals 

By 

Marie-Jeanne Munteanu and Paul Piraino 

The purpose of this paper is to apply the method studied in a previous 
paper (Munteanu [1]) for the three and four microwave channels which exist 
on TIROS-N. The context is exactly the same as in that paper, namely, the 
data used and the definition of the control experiment are identical. 

In the earlier paper we used 8 channels, 3 MSU and 5 IR. The reason for 
using only microwave channels in this paper is that the infrared channels used in 
the previous paper are affected by clouds in the case of real data. That 
study only dealt with simulated data, an idealized case. 

The results of the previous paper will have to be changed in a major way 
in order to be used on real data. We have applied the method of fuzzy 
partitioned clustering to predict tropopause height only using microwave 
information with an eye towards using it on real data under cloudy 
conditions. After this, in the second stage, we have developed stratified by 
tropopause regression temperature retrievals using only the three or four 
microwave channels for each 40 mb range. 

The data used are the 1575 radiosondes collected by N. Phyllips and their 
corresponding simulated microwave brightness temperatures. 

The first step in the experiment is the fuzzy partitioned clustering of 
the microwave brightness temperatures. This method is a combination of standard 
hard clustering and discriminant analysis. The new so called "fuzzy 
partitioned clustering" uses all the generated probabilities of membership of 
each pattern vector in any of the given clusters. These probabilities are 
generated by discriminant analysis in order to locate the correct cluster. 

The ultimate goal of standard discriminant analysis is to provide the 
unique (correct) cluster to which the pattern vector belongs. That means it 
takes into account only the maximum of all the generated probabilities and 
discards the rest of them. 

In our method we use all the probabilities and weight the regressions 
generated within each cluster. These regression formulas predict the 
tropopause height from the microwave brightness temperatures. With the 
standard method of hard clustering and discriminant analysis the overall 
r.m.s. error for prediction of the tropopause height is 42 mb for the summer 
and 49 mb for the winter. Using fuzzy partitioned clustering the overall 
r.m.s. error for prediction is 35 mb for the summer and 45 mb for the winter. 

In the second step of the experiment we have stratified the microwave 
regression temperature retrievals by tropopause height every 40 mb. The 
control experiment is defined as in the prior paper, namely, the data is 
stratified by land/ocean, summer/winter , and latitude bands. 
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A table with total r.m.s. errors for the stratified microwave temperature 
retrievals and the control experiment is displayed for three and four microwave 
channels. Significant improvement is obtained at most of the mandatory levels 
where we have satellite information using our method. 


Reference: 


1. Munteanu M. J., 1984: "The prediction of tropopause height from clusters of 
brightness temperatures and its applications in stratified regression 
temperature retrievals using microwave and infrared satellite measure- 
ments", (Proceedings of the Workshop on Advances in Remote Sensing 
Retrieval Methods , Williamsburg Nov.), submitted for publication at JCAM. 


Comparison of R.M.S Retrieved Temperature Errors Stratified by Tropopause 


1 Pressure 

Control 
Experiment 
4 MSU 

Stratified 
Regression 
4 MSU (40 mb) 

Control 
Experiment 
3 MSU 

Stratified 
Regression 
3 MSU (40 mb) 

30 mb 

2.12 

1.74 

2.11 

1.78 

50 mb 

2.21 

1.55 

2.23 

1.63 

70 mb 

2.02 

1.51 

2.09 

1.55 

100 mb 

1.94 

1.65 

1.92 

1.69 

150 mb 

2.25 

1.60 

2.41 

1.73 

200 mb 

3.00 

2.74 

3.13 

2.77 

250 mb 

2.69 

2.52 

2.77 

2.60 

300 mb 

2.69 

2.21 

2.76 

2.26 

400 mb 

2.91 

2.25 

2.90 

2.25 

500 mb 

2.75 

2.43 

2.73 

2.43 

700 mb 

2.31 

2.10 

2.51 

2.51 

850 mb 

3.28 

3.14 

3.33 

3.14 

1000 mb 

2.84 

2.92 

3.73 

3.33 

Total 

2.57 

2.19 

2.70 

2.35 
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Applications of Fuzzy Set Theory to Satellite Soundings 
Marie - Jeanne Munteanu 


Cao Hongxing and Chen Guofan have used fuzzy partitions based on 
imprecise relations to partition weather processes [ 1 ] . These authors have 
applied the principles of fuzzy similiarity to the partition of daily 
circulation processes during winter-spring of 1972. 

The objective of their study has been to group all circulations into 
different patterns that possess clear features in order to provide a basis 
for weather forecasting based on principles of fuzzy relations and cluster 
analysis. 

The purpose of this paper is to propose, first, the introduction of an 
appropriate fuzzy setting for satellite soundings and then, with this setting 
in mind, to apply it to clustering methods via unimodal fuzzy sets in the 
future. Methods of hard clustering analysis and fuzzy partitioned clustering 
have been applied by the author on simulated data with very encouraging 
results [3]. 

The new methodology proposed here is meant to apply to real data. 

We will start with a discussion of the clustering technique proposed. The 
notion of a unimodal fuzzy set has been chosen to represent the partition of a 
data set for two reasons. First, it is capable of detecting all the locations 
in the vector space where there exist highly concentrated clusters of points, 
since these will appear as modes according to some measure of "cohesiveness. " 
Second, the notion is general enough to represent clusters that exhibit quite 
general distributions of points. 

The technique is optimal in the sense that it detects all of the existing 
unimodal fuzzy sets and realizes the maximum separation among them. It is 
economical in memory space and computational time requirements and also detects 
groups that are fairly generally distributed in the feature space. The 
algorithm is a systematic procedure instead of iterative in nature. Most 
of the other algorithms are not capable of detecting categories with com- 
plicated distributions in the feature space and a great many are not appli- 
cable to large data sets. 

An important distinction between the clustering unimodal fuzzy algorithm 
and other clustering algorithms is that the latter use a distance measure 
as the only means of clustering. In [2] Gitman and Levine have 
introduced another "dimension", the dimension of the order of importance 
of every point, as an aid in the clustering process. This is accomplished 
by associating with every point in the set a grade of membership or charac- 
teristic value. Thus, the order of the points according to their grade of 
membership as well as their order according to distance are used in the 
technique. The order of importance is introduced via the grade of member- 
ship notion which, in turn, is introduced via the fuzzy set setting. 

This latter partitions the sample from a multimodal fuzzy set into unimodal 
fuzzy sets. This is a very important step since any inversion method will 
perform considerably better once we divide the multimodal data sets into 
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ones that are uniraodal in nature. The motivation behind this new methodo- 
logy is the fact that our satellite data as well as our ultimate goal, the 
set of temperature profiles, are multimodal in nature. 

Now we will present the fuzzy sets setting which we propose for the 
satellite soundings, a setting which has merits by itself and which will 
also provide the context for the unimodal fuzzy clustering method proposed. 

Fuzzy Sets 

The theory of fuzzy sets has been invented with the idea of introducing 
a multivalued logic which could possibly simulate some processes of the 
human brain. It is well known that the human brain performs much better 
than any computer in the world in the field of pattern recognition (and 
many other subtle qualitative decision-making processes). 

Experts in artificial intelligence have only begun to use fuzzy sets. 

This context seems appropriate for satellite sounding because of the nature 
of the information being associated with a specific weighting function. 

Consider the weighting functions as membership functions and form 
pairs of the soundings and their corresponding membership functions. There 
are two very interesting cases. First, take the satellite information 
at the peak of the weighting functions or as representing the information 
at the peak of the weighting fuctions. In terms of fuzzy sets this repre- 
sents reducing or approximating the weighting function by one point, namely 
its peak. The second important case is the mean layer case. We take 
slices of the weighting functions contained within the layer. There it is 
possible to define how to extract the optimal information within a layer 
as well as intermediary steps of approximations. In this way we can define 
maximum information within a layer and eliminate redundancy by the nature of 
operations in the fuzzy sets context. 

This new context has the potential of providing solutions for the optimal 
distribution of channels in order to solve important structure of the tempera- 
ture profile and optimal ways of scanning given certain constraints. For this, 
methods of optimization theory (among them dynamic programming) should be used. 

For the first case, when information from a given channel is considered 
at one point, classification of sounding data has been proposed with hard clus- 
tering and fuzzy partitioned clustering on simulated data. These techniques 
have been the object of two papers published by the author in the Williamsburg 
1984 proceedings. 


Further generalizations 

If we consider the N-tuple element containing the satellite measurements 
from one given channel and its varying membership function (with the angle, 
the angle time, and any other important quantities) in the general context 
of distributive lattices, then we can define different optimal methods of 


34 


extracting information from our satellite data. 

For these more complicated assumptions we have to generalize fuzzy 
sets theory to fit our needs while keeping the same basic assumptions. 
Specifically, we have to consider fuzzy sets in the framework of distributive 
lattices. Once the optimal information is derived from satellite data, 
different methods of inversion could be applied. We also can have an 
estimate of the error by taking different degrees of approximation of the 
given satellite data. 

Another important point is that combining data optimally could lead to 
new inversion methods. For example, if we take besides the satellite 
data, its membership function as a function of angle and the angle itself we can 
conceive an arrangement of data containing optimal information leading 
naturally to the use of the computerized tomography method (used by H. Flemming). 

One new important approach would be to arrange data optimally for 
different two dimensional retrievals. We may take information at grid points 
or along a satellite track. 

In more dimensions we may have a two dimensional field such as total ozone 
which, for example, we want to integrate in the general scheme as a constraint. 
Taking such surface constraints would represent an optimal approximation in 
three dimensions while taking information at a point would represent a minimal 
approximation. Taking data along satellite tracks or other curves in 
space would be an intermediate approximation for a three-dimensional solution. 

In the last part of the abstract a few definitions of fuzzy sets theory 
are presented. 

The theory of fuzzy sets 

By fuzziness we mean a type of imprecision which is associated with 
indistinct sets, that is, classes in which there are no sharp boundaries and 
no sharp distinction between membership and nonmerabership • 

Definition ; Let X = {x} denote a collection of objects (points) x. Then a 
fuzzy set A in X is a set of ordered pairs 


A = {(x,pa (x)): x e X} 

where y^(x) i- s termed the grade of membership of x in A and y^: X + M is a 
function from X to a space M called the membership space. 

The choice of M 

1. We can assume H as the interval [0,1] with 0 and 1 representing respectively 
the lowest and the highest grades of membership. 

2. More generally, M can be a partially ordered set or more particularly a 
lattice (which has a maximal and minimal element). 
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Thus our basic assumption is that a fuzzy set, despite the unsharpness 
of its boundaries, can be defined precisely by associating with each object x 
a number between 0 and 1 which represents its grade of membership in A. 

Definitions 


Normality : A fuzzy set A is normal iff 

Su Px PA<x) = 1. 

Equality : Two fuzzy sets A and B are equal iff y A = y B that is 


y A (x) = ygCx) 


for all x e X. 

Unntainment : A fuzzy set A is contained in or is a subset of a fuzzy set B, 

A c B iff 


PA < MB* 


Complementation : A' is said to be the complement of A iff 


PA' = 1 ~ PA* 

Intersection : AOB is said to be the intersection of A and B iff 

MAftB (x) = Min(y A (x) , y B (x)) 

for all x e X. 

Union: A (J B is said to be the union of A and B iff 

PAUB (x) = Max (p A (x), p B (x) ) 

for all x e X. 

Union and intersection are associative, they distribute over each other and 
satisfy De Morgan's laws. It is easy to verify that 

MA' U B' = P(AOB) ' 

MA’O B' = P(AUB) ' 
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Algebraic Product: 


AB is said to be the algebraic product of A and B iff 

hab( x ) = v*a( x ) * ub( x ) 

for all x s X. 

Algebraic Sum : 

A + B is said to be the algebraic sum of A and B iff 

- J A+B( X ) = UA< X )+ ! j b( x ) “ Pa( x ) * Pb( x ) 

for all x e X. 


The algebraic sum and product are associative, but do not distribute over 
each other. They satisfy a form of De Morgan's laws. It is easy to verify 

(A + B)' = (A'B*)'. 

The fuzzy sets as described here as well as the unimodal fuzzy clustering 
are introduced with an eye towards developing an expert system for the 
processing and interpretation of real satellite data. 
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Applications of Some Artificial Intelligence Methods to Satellite Soundings 
Marie- Jeanne Munteanu and Oleg Jakubowicz 

In this paper we continue the research started in a previous paper 
Munteanu et. al_.[l] The new element is the use of fuzzy partitioned clustering 
instead of hard clustering for temperature profiles. 

We use exactly the same setting, namely, 1575 radiosondes collected by 
N. Phyllips and the corresponding simulated satellite measurements. The control 
experiment is formulated in exactly the same way, namely, the usual concept of 
regression used by NMC based on stratification by land, ocean, summer, winter 
and latitude bands. 

In the first paper we used hard clustering of temperature profiles and 
discriminant analysis to locate the correct cluster from satellite measurements. 

In the second step we developed regression temperature retrievals within 
each cluster and compared the total r.m.s. errors against the control experiment. 
Excellent results have been obtained, namely, improvements of more than 1° 
appeared at all of the mandatory levels where we had satellite information. 

In this paper our aim was to refine this method using the probabilities 
of membership of each pattern vector in each of the clusters derived with 
discriminant analysis • 

In hard clustering one takes the maximum probability and considers the 
corresponding cluster as the correct cluster discarding the rest of the 
probabilities. In fuzzy partitioned clustering we keep these probabilities and 
the final regression retrieval is a weighted regression retrieval of several 
clusters. We have used this method in the clustering of brightness temperatures 
where the purpose was to predict tropopause height. (See [2]). 

This is only the first step towards a fuzzy clustering where we minimize 
more qualitative functionals than just the euclidean distance. A further 
refinement in this work is the division of temperature profiles into three major 
regions for classification purposes. 

The results are summarized in the tables where we display total r.m.s. errors. 
It seems that an approach based on fuzzy logic which is intimately related 
to artificial intelligence methods holds promise. Larger improvements than in 
the case of hard clustering are obtained at most of the levels. 
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Comparison of RMS Temperature Retrievals (Summer) 5IR & 3 MSU 



Figure 1 


Comparison of RMS Temperature Retrievals (Winter) SIR & 3 MSU 


Pressure 

Control 

Experiment 

Clusters of 
Temperature 

Fuzzy Partitioned 
Clusters of 
Temperatures 

30 mb 

4.71 

4.19 

4.09 

50 mb 

2.88 

2.43 

2.35 

70 mb 

1.89 

1.74 

1.70 

100 mb 

2.31 

1.94 

1.89 

150 mb 

2.15 

1.97 

1.88 

200 mb 

2.36 

2.16 

1.70 

250 mb 

2.43 

1.86 

1.80 

300 mb 

2.53 

1.81 

1.69 

400 mb 

2.40 

1.73 

1.51 

500 mb 

2.62 

1.77 

1.63 

700 mb 

2.57 

1.86 

1.69 

850 mb 

2.75 

2.18 

2.01 

1000 mb 

2.72 

2.30 

2.10 


Figure 2 
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C. ANALYSIS AND FORECAST MODEL DEVELOPMENT 



A Factored Implicit Scheme for Numerical Weather 
Prediction with Small Factorization Error 

J. M. Augenbaum, S. E. Cohn 1 , and D. Marchesin 2 


The time step used in current models for integrating the global baroclinic 
primitive equations is limited by the Courant— Friedrich— Lewy (CFL) stability 
criterion. In schemes based on explicit time discretization, for example, 
whether finite-difference, pseudospectral or finite element in space, the time 
step cannot exceed the quotient of the smallest spatial interval and the maximum 
phase speed. This is a severe restriction first of all because the fastest 
waves propagate roughly ten times faster then the waves which carry most of 
the atmosphere's synoptic scale energy. Second, on a uniform latitude-longitude 
grid the mesh spacing becomes very small near the poles. As a consequence the 
time discretization error of explicit schemes is much smaller than the spatial 
discretization error, especially in middle and low latitudes. 

By using a scheme with more desirable stability properties, it is 
possible to increase the time step without effecting overall forecast accuracy, 
thereby increasing computational efficiency. Recently, Cohn et al. (1985) have 
introduced a fully implicit scheme for the barotropic primitive equations based 
on the work of Beam and Warming (hereafter referred to as B-W). The B-W scheme 
is unconditionally linearly stable and allows use of a large time step, whose 
size is determined solely by accuracy requirements. The method is based on a 
spatial factorization of the implicit operator into one— dimensional operators, 
and yields second order accuracy in time and fourth order accuracy in space. 

The scheme is efficient because it only requires the solution of one-dimensional 
block-t rid i agonal linear systems. In an effort to analyze the accuracy of the B-W 
scheme, for meteorological applications, we have carried out an analysis of the 
factorization error for a linear shallow water system on a tangent plane. By 
discretizing the equations only in the time domain, while assuming a spatial 
variation of the form W™ = W™ e*' 1 ** + ly) we can derive analytic expressions 
for the amplification matrix, F(k,l), in terms of the wave number (k,l) i.e. 

l) n+1 = F(k,l) 

The eigenvalues of F, for a given set of wave numbers, can then be used to find 
amplitudes and periods of the various wave components of the solution. This 
method can also be used to find the periods and amplitudes of the various modes 
for the unfactored Crank-Nicolson scheme (referred to as C-N) and, with a slight 
modification, the continuous (true) solution (referred to as CONT). 

Numerical results (presented in table 1) show that, for large time steps, 
the factorization error can be significant, even for the slowly propagating 
Rossby modes. We have formulated a new scheme based on a more accurate 
factorization of the equations. By grouping separately the terms of the 
equations which give rise to the fast and slow motion, the equations are factored 

1 Courant Institute of Mathematical Sciences, New York University 
New York, N. Y. 10012 

2 Department of Mathematics, Pontificia Universidade Catolica do Rio de 
Janiero, Rio de Janeiro-RJ, CEP22453, Brazil 
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in a more accurate way. The fast-slow factorization completely eliminates the 
factorization error. If each of the fast and slow factors are now factored 
again, this time according to spatial components, the resulting scheme only 
involves the solution of one-dimensional linear systems, and is computationally 
efficient. We have shown that the factorization error for the slow mode 
component of the solution is negligible for this new scheme (which we refer to 
as F-S). 

In table 1 we present a comparison of several implicit schemes (including 
the semi-implicit scheme, which we refer to as SIMP) for a Rossby wave with 
wavenumber 16,8 in x,y directions. 

The results of table 1 show that the fast-slow scheme is accurate, even 
for large time steps, for middle and low latitudes. 

We are currently conducting numerical experiments to compare the various 
schemes for the full nonlinear barotropic shallow water equations. 

Table 1 

Period (in hours) for the Rossby wave component of a 16 x 8 wave with 
At = 30, 60, 120 min. at lat 0 = 0°, 45°. 



Lat = 0° 

Lat = 45° | 

At (min) 
scheme 

30 

60 

120 

30 

60 

120 

CONT 

27.82 

27.82 

27.82 

20.90 

20.90 

20.90 

C-N 

27.85 

27.94 

28.30 

20.94 

21.08 

21.59 

B-W 

30.70 

39.35 

74.01 

23.60 

31.73 

64.43 

F-S 

27.84 

27.89 

28.07 

20.93 

21.01 

21.33 

SIMP 

ro 

■ 

O' 

27.58 

26.81 

20.81 

20.53 

19.25 
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Comparison of Optimum Interpolation and Cressman Analyses 
W. E. Baker, S. C. Bloom, and M. S. Nestler 


The objective of this investigation is to develop a state-of-the art 
optimum interpolation (0/1) objective analysis procedure for use in numerical 
weather prediction studies. A three-dimensional multivariate 0/1 analysis 
scheme has been developed for the purpose. Some characteristics of the GLAS 
0/1 compared with those of the NMC and ECMWF systems are summarized in Table 
1. Some recent enhancements of the GLAS scheme include a univariate analysis 
of water vapor mixing ratio, a geographically dependent model prediction error 
correlation function and a multivariate oceanic surface analysis (Bloom et al., 
1984). 


Figure 1 compares the 24 h and 48 h forecasts from 0000 GMT 13 January 1979 
initial conditions provided by the 0/1 analysis scheme with forecasts from a 
successive correction method (SCM) of analysis developed earlier. As may be 
seen in Fig. 1, the forecast of the intensity and position of the cyclone which 
moves northeastward across the United States as well as the ridging which follows, 
is more accurate with the 0/1 than with the SCM. 

Future work includes completing vectorization of the 0/1 on the Cyber 205, 
a detailed comparison with the SCM, real and simulated data impact studies and 
experiments with normal mode initialization. 
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TABLE 1. Characteristics of the GLAS 01 compared with those of other systems. 


SALIENT 

FEATUKES 

GLAS 

NMC 

ECMWF 

APPROACH 

3-D; MULTIVARIATE IN Z,U,V 
UNIVARIATE IN W 

3-D; MULTIVARIATE IN Z,U,V 
UNIVARIATE IN RH 


RESOLUTION 

18 LEVELS; 4° x 5° (2-U) 
12 LEVELS; 4° x 5° (3-D) 

12 LEVELS; 3-75° UPPER AIR, 
2-5° SURFACE 

15 LEVELS; 1.875° 

AUTOCORRELATION 

FUNCTION 

DAMPED COSINE FOR Z 
DAMPED EXPONENTIAL FOR W 

6AUSSI AN 

BESSEL 

CROSS- 

CORRELATION 

MODEL 

GEOSTROPHIC (SCALED TO 0-0 
AT EQUATOR) 

SAME 

SAME 

ASSUMPTION OF 
ISOTROPY 

Z (YES); U,V (NO) 

SAME 

SAME 

GEOGRAPHICALLY 
DEPENDENT v 

YES 

NO 

1 

YES (LATITUDINALLY 
DEPENDENT) 

MULTIVARIATE 
SURFACE ANALYSIS 

YES (EKMAN BALANCE) 

YES (GEOSTROPHIC) 

YES (GEOSTROPHIC AT 
1000 MB) 




















SCM 


Observed-14 January 


24 h 

Fcst 



SCM 


0/1 


Observed-15 January 



Fig. 1. 01 and SCM forecasts from 0000 GMT 13 January 1979 






Derivation of Model Topography 
R. C. Balgovind 


In the GLA Fourth-Order model it is necessary to represent the topography. 

In doing so one has to address the problem of the representation of the topography 
at grid points. Topographies we prepared were for the 4° x 5° and 2° x 2.5° 
forecast models. The tabulated source data came from the Navy Fleet Numerical 
Oceanography Center at Monterey. This source has the global elevation given 
at every 10 minutes in latitude and longitude. The accompanying figures show 
six examples of different topographys [01, SI, S23, S3, Q3 and Tl] for the 
4° x 5° model derived from the source data. 

We first attempted to derive an envelope topography. An example of this 
is the Tl_ obtained by taking local mean plus one standard deviation at each 
grid point and sigma filtering it. We found this method to be greatly influenced 
by large standard deviations at steep mountains. 

The G1 topography is just the local mean. The S]^ is obtained by Sigma* 
filtering in both latitude and longitude the mean 01 . S2 is when the operation 
is applied twice and jS3^ thrice. Q3 is the sigma filtered local mean of the 
upper third quantile of the source data - that is the mean of the one third of 
the largest source data. 

Presently experiments are underway to determine which is the best topography 
for the Fourth-Order model. Preliminary results indicate a close contest between 
Q3 and SI. 


* For Sigma Filter see paper by that name in this review. 
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Sigma Filter 


R. C. Balgovind 


In deriving the topography for the GLA Fourth-Order model we needed to 
smooth the topography. This is in order to remove the Gibbs phenomenon. The 
Gibbs phenomenon occurs whenever we truncate a Fourier Series . 

Sigma filtering is smoothing using Sigma factors. Cornelius Lanczos (1966) 
introduced the Sigma factors for the purpose of reducing the Gibbs phenomenon. 
The idea is that for a given Fourier Series 

N 

ffl(t) = a 0 /2 + £ [afc cos kt + b^ sin kt] 

k=l 


we apply the averaging operator 


_N 

2n 


t + t/N 

/ f N (x)dx 

t - t/N 


The result of this is a smooth field 

a N 

f n (t) = — + ]> [a^ coskt + b^sinkt] 

2 k=l 



That is the smooth, Fourier series is nothing but the orginal Fourier series 
with its coefficients multiplied by corresponding sigma factors. 

This operator can be applied many times to obtain high order sigma filtered 
field and is easily applicable using FFT. We found this filter beneficial in 
deriving the topography. 
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Recent Developments in Nonlinear Normal Mode Initialization 


S. Bloom 


In order to assess the importance of a balanced initial condition upon GLAS 
GCM forecasts and assimilation cycles, an effort to combine previous work on 
normal mode initilization at GLA (Bloom, 1982, 1983; Navon et al . , 1983) is 
currently underway with the goal of developing an initialization process for 
the production version of the GLAS 4th order GCM. The major aspects of this 
work fall into two parts: vectorization of the linear projector code (Bloom, 

1983), and the insertion of the mode projector and Machenhauer iteration algorithm 
into the full GLAS GCM. 

Memory and paging constraints place restrictions on the number of horizontal 
modes stored for initialization purposes, as well as on the manner in which they 
are stored. Only the first five vertical structures of the gravity modes are 
used; these account for only 40% of the total set of the model’s normal modes. 

One source of difficulty encountered in moving the projector code to the CYBER 
205 lay in the use of different FFT packages on the two machines (NCAR on the 
Amdahl, Temperton on the Cyber). Differing phase and normalization conventions 
provided many elusive coding errors. 

For the preliminary studies, a Machenhauer nonlinear normal mode initializa- 
tion technique (Williamson and Temperton, 1981) is used. This method entails 
the insertion of a modified version of the mode projector into the full GCM, 
and the modification of the GCM to allow for iterative calls to the projector. 

The main difficulties encountered with inserting the projector code have been 
in obtaining the correct communication among the various model and projector 
commons, and in performing the update to the correct variables within the 
projector. The GCM is integrated for one time step in order to estimate the 
time tendencies of the ’fast’ gravity mode coefficients; however proper care 
must be taken to save the original set of fields since the ’slow’ coefficients 
are never computed. If the original fields are not saved, the iterative process 
diverges. This result is due to the ’slow’ coefficients evolving in time away 
from the adjusted ’fast’ coefficients; the resulting unbalanced fields quickly 
develop unrealistic winds. The GLAS GCM also has a special treatment of the 
fields at the poles (using subroutines POLINP and P0L0UT) ; if this is ignored 
during the iterations, wildly unrealistic winds and pressures ensue. Finally, 
if a timestep is chosen too small (for the model integration to estimate time 
tendencies), then truncation error will cause the iterations to diverge. 

The initialization process is still undergoing testing. Preliminary 
results do show the Machenhauer process to be convergent over 4 iterations if 
only the external vertical structure is initialized. However, tests initializing 
the first two vertical structures diverge after the second iteration. Currently, 
a shallow water version of the GCM (the same one as used in Navon et al . , 1983) 
is being used to test various parts of the initialization algorithm. 
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Comparison of Two Orographical Data Sets for the GLAS Fourth Order GCM 


H. M. Helfand, R. Balgovind, R. Dlouhy, J. Pfaendtner 


The development of a 2° latitude by 2.5° longitude version of the GLAS 
Fourth Order General Circulation Model has necessitated the specification of a 
correspondingly fine horizontal resolution orographic data set to act as a lower 
boundary condition for the GCM. Elsewhere in this volume, R. Balgovind has 
described the development of two such fine resolution orographic data sets 
which he has labeled the SI and the Q3 topographies. These differ from one 
another primarily in their respective height specifications over rugged terrain. 

The SI topography attempts to relate the mean areally-averged height of such 
terrain while the Q3 topography emphasizes the effect of the flow of the tallest 
peaks in such regions by specifying an exaggerated "significant height" as 
described by Balgovind. 

In order to assess the relative merits of these two orographic data sets, 
we have initiated a series of tests to evaluate the comparative weather forecasting 
and climate simulation skills of the GLAS Fourth Order GCM when the SI and Q3 
topographies are applied to the lower boundary. These tests were carried out at 
both 2° x 2.5° and 4° x 5° horizontal resolutions. 

Three sets of Northern Hemispheric winter forecasts at the 2° x 2.5° 
resolution indicate that the Q3 topography gives slightly better results over 
the Northern Hemisphere than does the SI topography. The improvement did not 
occur until after 4 days for Skill Score 1 and after 2 to 3 days for root mean 
square error. Fig. 1 illustrates a typical example of this improvement. There 
were no significant trends for forecasts over the Southern Hemisphere on tropical 
regions. No forecasts were available at this resolution for the summer season. 

The results of the 4° x 5° forecasts were mixed. For the first 3 to 4 days, 
sea level pressure forecasts for the Northern Hemisphere were better with the 
Q3 topography for a period of time. After this time, there were no consistent 
trends. 

In the Southern Hemisphere, forecasts with the Q3 topography were slightly 
better than those with SI topography but not until after 3 to 6 days. In the 
tropics, a June forecast was substantially better with the SI topography while 
the (Northern Hemisphere) winter forecasts favored the Q3 topography at sea 
level and were mixed at 500 mb. 

It is apparent that more forecasts will be necesary at both horizontal 
resolutions before any definitive conclusions can be made regarding the choice 
of an optimal orographic data set. 
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Dotted = 2° x 2*5° with SI topography 
Ticked = 2° x 2.5° with Q3 topography 
Dashed = 4° x 5° with SI topography 
Heavy = 4° x 5° with Q3 topography 
Solid = 4° x 5° with old topography 
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Fig* 1 Forecast skill for the Northern Hemisphere 
for a forecast beginning on Jan. 5, 1979. 

The top panels represent sea level pressure 
scores, while the bottom panels represent 
500 mb heights. Skill score 1 is shown 
in the left panels, and root mean square 
error on the right. 










Development and Testing of the Variable 
Vertical Resolution Fourth Order GCM 


H. M. Helfand, R. Dlouhy, J. Pfaendtner, L. L. Takacs 


The vertical coordinate of the Fourth Order Model has been generalized so 
that the model can now run with an arbitrary number of vertical layers and so 
that the thicknesses of these layers can be arbitrarily specified (in the o- 
coordinate). This Variable Vertical Resolution (VVR) version of the Fourth 
Order Model will soon replace the current production model although it will 
probably still be run with 9 equally spaced layers when used for production 
purposes. 

In addition to the generalization of vertical interpolation and quadrature 
calculation and of vertical array dimensioning, the VVR model contains the 
following changes: 

i) The strapping scheme for cumulus convection has been generalized from 
the current 2-2-2 scheme to a three-layer scheme whose strapping para- 
meters can be arbitrarily specified through namelist parameters. 

ii) The redistribution of convective heating is now specified by a con- 
serving quadratic interpolation scheme and the scheme for the redistri- 
bution of convection moistening has been generalized. 

iii) A relative humidity criterion for the onset of cumulus convection 
(Helfand, 1981) has been specified as a namelist option. 

iv) Subroutine C0MP3 has been fully modularized for ease of testing and 
development of new physical parameterization schemes. 

v) The computation of large-scale precipitation has been improved by 

limiting the evaporation of falling rainwater to the amount which can 
be held by the receiving atmospheric layer as is computed by linear- 
ization about the state of that layer before such evaporation occurs. 

vi) Optical thicknesses of clouds are computed to reflect changes in the 
vertical grid and in the strapping parameters for convective clouds. 

vii) Planetary Boundary Layer (PBL) parameters are computed in the basis 
of the depth of the lowest model layer. 

Presently physical parameterization schemes are being developed so that the 
vertical structure of the PBL can be explicitly resolved in the VVR model but this 
work is currently to be regarded only as experimental. 

To assess the skill of the WR model, it has been run with 9 equally spaced 
layers and compared with the current production model. In two Northern Hemis- 
pheric winter cases and one summer case, the two models were virtually identical 
in forecast skill for 6 to 7 days. After that the VVR model was slightly better 
in the winter cases and the production model was slightly better in the summer 
case. The only exception to this was that after 2 days the production model 
gave slightly more skillful 500 mb forecasts in the tropics for the summer case. 
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Surprisingly, there were substantial differences in the climatologies 
produced by the two models. The mean January 500 mb height fields for simula- 
tions initialized at OZ, 15 December 1978, for example (see figure 1), were 
markedly different in the location and shape of a North Atlantic ridge, in the 
location of a trough over central Europe, in the shape and westward extent of 
the Siberian low, and in the zonality of the Southern Hemispheric and mid-lati- 
tude westerlies. 

The simulation by the VVR model seems to be more in accord with mean 
observations for January 1979 which suggests that the VVR model might give better 
extended range weather forecasts. The simulation by the production model, on 
the other hand, appears to be in better agreement with the climatological mean 
with the exception of European trough. Obviously, the comparison of the two 
models warrants further study. 
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Mean geopotential height field for the period 1 January - 1 February for a) an integration with 
the VVR model begun at 0Z 15 December, 1978, b) a similar integration with the production model, c) 
observations for the month of January, 1979, and d) climatology for the month of January. 









Parameterization of Surface Fluxes 
in the VVR Fourth Order GCM 

H. M. Helfand 


The Variable Vertical Resolution (VVR) option of the GLAS Fourth Order GCM 
(described elsewhere in this volume as well as in the 1980/81 Research Review) 
allows one to enhance the vertical resolution of the region of the atmosphere 
adjacent to the earth's surface* This, in turn, makes it possible to compute 
turbulent suface fluxes of heat, momentum and moisture directly from the prog- 
nosticated properties of the lowest model layer by use of the Monin-Obukhov 
surface layer similarity theory. However, one must be extremely cautious in 
doing so, for the similarity theory applies formally only to the "constant 
flux" surface layer which is but a few tens of meters deep. It is not practically 
feasible to work with a lowest GCM layer thin enough to satisfy these formal 
constraints. 

It is fortuitous that, at least under conditions of neutral stratification, 
the similarity theory can be extended beyond its formal limits. Panofsky (1973) 
points out that, because of the cancellation of opposing effects (decreasing 
stress and less-than-linearly increasing mixing length), the theory still gives 
reasonable looking results, when the distance from the ground becomes as large 
as 150 m. ° 


It is not infeasible to run the VVR model with a lowest layer thickness 
on the order of 300 m (the center of the layer would then be 150 m above the 
earth's surface), and so if one can prescribe similarity functions <|> ra (s) and 
<f>h(0 to adequately describe the entire "extended surface layer", the problem 
of surface flux parameterization is solved. 


For small values of |?| , the similarity functions should approach the 
behavior of the similarity functions of conventional surface layers such as 
those of Businger et al. (1971): 


= 

i - v> _1/ * 

C < o 



1 + 

c 1 0 


( 

r a - v hC r 1/2 

5 < 0 

(1) 

<f>h(S) = *v 



l 

. 1 + «hC 

5 < 0 



As |?| becomes large, <j> m and should exhibit the appropriate asymptotic 
behavior. For the unstable case C < 0, this behavior is 


^m,h^) 


.-1/3 


( 2 ) 


For the stable case £ > 0, Clarke's (1970) mean stable case for the Wangara 
Experiment suggests the asymptotic form 
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-1 


( 3 ) 


To interpolate between the limits (1) and (2), the similarity functions 
for the unstable case 5 < 0 have been taken as the solutions to 


‘t’m ~ - 1 - 

‘t’h 2 “ ^h^^h 3 = 1 * 


For the stable case C > 0, the similarity fuctions are 

^m,h^^ 


1 + a m,h C 1 


1 + 3m , \fi (1 + <*m,hCl) 


(A) 


(5) 


where 


* 

Cl = min (Cl, C m ,h) • 


Thus (5) interpolates between (1) and (3). 


For the stable case. 


W c) 


can be integrated directly, while for the 


unstable case, the integrals can be integrated in terms of the solutions of (4) 
which can be solved, in turn, through mathematical iteration. The momentum and 
heat transfer coefficients so obtained are shown in Fig. 1. Note that they 
approach the solutions of Businger et al . (1971) for small I Rib I • ^ 

Rib * ~°° one has 


C H - ^T ~ Rl b^ 2 


which gives the proper asymptotic limit 

W0 S = C R U AO - A0 3 / 2 . 


For Rib -*■ 00 the surface fluxes remain positive which correctly reflects 
the facts that the Richardson number remains subcritical in an increasingly 
narrower neighborhood of the earth's surface. 
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Fig. 1 Momentum and heat transport coefficients Cy and Op versus bulk Richardson 
number Riy. Solid lines represent present model while dashed lines 
represent the model of Businger et al . (1971). Graphs are given for 
h/ Z q equal to 10, 10 4 and 10 5 where h is the depth of the surface 
layer and zq is the surface roughness height, described elsewhere in 
this volume. The inset shows the behavior of the curves on an expanded 
Richardson number scale about the neutral case Riy = 0. 
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Specification of Surface Roughness 
Over Oceans in the GLAS Fourth Order GCM 


H. M. He If and 


The surface roughness height zq is a parameter which measures the 
effect of surface irregularities on the mean profiles which occur in the 
atmospheric surface layer. Over land surfaces, zq is a function of soil 
composition, plant canopy type and distribution and other factors such as 
buildings and structures. Over open water, zq depends only on the 
distribution of waves, wavelets and other surface irregularities. These, 
in turn, are thought to be determined by the mean wind within the surface 
layer. 


Large and Pond (1981) have derived an empirical relationship 
between the mean wind speed Ujo at the ten meter level over open water 
and the neutral drag coefficient Cpfj Their relationship covers the range 

Ujq _> 2 ms - ^ . 

Earlier, Kondo (1975) had computed a relationship between Uiq 
and Cj)N on the basis of laboratory observations which extended down to 
the range of infinitely small wind speeds. Fig. 1 is an interpretation 
between the relationships of Large and Pond for large wind speed and of 
Kondo for small wind speed. 

One can convert the relationship of Fig. 1 to the relationship 
between the wind stress velocity 


u*= 



U 


10 


and the surface roughness length zq. This has been done for the Monin - 
Obukhov surface flux scheme in the GLAS Model so that zq can be expressed 
as 


= — + A 0 + Ao u. + A/, u. + Ac u 


u 


Where the A^ are defined separately over three separate domains of u . 
The function (1) is shown in Fig. 2. * 
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NEUTRAL DRAG C0EFFICIENT — CDN 



SURFACE R0UGHNESS LENGTH (CM) 



Fig. 1. Neutral drag coefficient 

expressed as a function of 10m 
wind speed U^q given in ms 
This function interpolates 
between that of Large and Pond 
(1981) for strong winds and 
that of Kondo (1975) for weak 
winds. 


Fig. 2. Surface roughness height zq ir 
cm expressed as a function of 
the surface stress velocity u* 
given in ms“l. This function j 
computed on the basis of the 
relationship shown in Fig. 1. 
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Development of an Optimum Interpolation Analysis Method for the CYBER 205 

M. Nestler, J. Woollen, and Y. Brin 


The GLA optimum interpolation (01) analysis method is a state-of-the-art 
technique to assimilate the diverse observational database obtained during FGGE, 
and thus create initial conditions for numerical forecasts. The 01 analyzes 
pressure, winds, and temperature at sea level, mixing ratio at six mandatory 
pressure levels up to 300 mb, and heights and winds at twelve levels up to 50 mb. 

Conversion to the CYBER 205 required a major re-write of the Amdahl 01 
code to take advantage of the CYBER vector processing capabilities. Structured 
programming methods were used to write the programs and this has resulted in a 
modular, understandable code. We note here that the conversion effort is 
continuing at the time of this writing. 

We have included some degree of flexibility in the vector code in anticipation 
of future experiments with the analysis. The options available include: 

1) variable horizontal resolution (4° x 5°, 2° x 2.5°) 

2) variable vertical resolution (12 or 18 levels) 

3) a two-dimensional or three-dimensional analysis 

4) a two-dimensional troposphere/stratosphere analysis to 0.4 mb on a 
4° x 5° x 18 level grid 

5) the use of significant level data 

6) the use of up to 20 observations to influence the analysis of a grid point 

7) the use of different forecast error correlation models 

Table 1 compares the Amdahl and CYBER CPU time required for one three- 
dimensional analysis on a 4° x 5° x 12 level grid. A maximum of twenty observations 
influences any grid point. Among the contributors to the increased speed of 
the CYBER code are a vectorized covariance-calculation routine, an extremely 
fast matrix equation solver, and an innovative data search and sort technique. 

The solver stores up to 400 matrix equations for simultaneous solution. During 
matrix decomposition, time is saved by computing the expensive reciprocal of 
the square root operation on only selected diagonal terms. Vectorized routines 
move the matrix elements into and out of the format needed for solution. To 
locate data for use in the analysis, all grid points in the vicinity of a data 
point are listed. After each data point has an associated list of gridpoints, 
the lists are sorted so that each gridpoint now has an associated list of data 
points. This search and sort routine has also contributed to a marked decrease 
in wall clock time required for the analysis. 
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Model Development Highlight for 1984: 
The GLA 4th Order GCM 

J. Pfaendtner 


A number of improvements have been made to the GLA Fourth Order GCM and 
its associated run procedures during the past year. The purpose of this note 
is to document the major changes which have been made, and to indicate the 
current direction of our model development efforts. 

• Model Resolution 

Model versions using the 4 degree latitude by 5 degree longitude grid 
resolution and 9 sigma layers are now completely core contained in the 1 M word 
CYBER memory. They need 175 seconds of real time to complete a standard 1 day 
forecast (Matsuno time scheme, 7 1/2 minute timestep, physics every 30 minutes, 
longwave radiation every 3 hours, history records written at 6 hour intervals). 
Forecasts are also routinely being done with a 2 degree latitude by 2 1/2 degree 
longitude, 9-layer version of the model. At this resolution the hydrodynamics 
is also core contained, but paging occurs in the model physics. A one day 
forecast requires 26 minutes of real time and 14 minutes of CPU time. A 15- 
layer version, which can use non— equidistant sigma levels, is also running in 
core at the 4 by 5 horizontal resolution. 

On April 1,1985 an additional 1 M words of memory are scheduled to be 
installed on the CYBER 205. With this additional space all model versions should 
run in core. It is projected that a planned 2 by 2 1/2 degree, 15-layer model 
will require 22 minutes of real time for one simulated day. The code management 
and run procedures have been modified so that any version of the model can be 
run at any desired resolution from a single set of resolution independent masters 
for that version. 

© Boundary Fields and Topography 

The model was modified to allow the boundary fields (sea surface temperature, 
albedo, ground wetness, sea ice extent) to be interpolated in time during a 
forecast or climate simulation. The fields can be taken either from climatology 
or from the actual observed anomalies during the period being simulated. This 
feature was used extensively for the recent study of the 1980 heat wave drought 
(Wolf son and Atlas). It is also being used for a study of the effect of local 
SST anomalies on the Ethiopian drought. A number of improved topography fields, 
including an envelope topography, have been prepared for both horizontal resolu- 
tions of the model (Balgovind). The effects of the altered topographic forcing 
on medium range forecasts and the model climatology were considerable (Helfand). 

© Post-processing of model results on the CYBER 

Changes were made to the model to enable it to output history fields on 
mandatory pressure levels as well as on sigma surfaces. This feature, which 
will appear in the next production version of the model, is important in two 
respects. First, it will provide an orderly sequence of pressure level history 
tapes for each experiment with data management handled automatically by the 
model merge job. Second, when running the forecast/analysis cycle, it will 
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enable the direct examination of the analysed fields. Until now, the GLAS 
analysis has been interpolated from the pressure levels on which it is done to 
the model sigma levels for inclusion in the sigma level history tapes. An 

the 6 analysis 11 fr ° m ^ Sigma levels back t0 pressure was required when examining 


^ Th pvfuD 6 u ati v n ° f verification statistics using the pressure level history 
on the CYBER has been implemented. Plans are being made to extensively expand 
the options available for examination of model results on the CYBER. Model 
energetics and time averages will be the first options to be added. 


° Improvements to model source code 


The model code for parameterizations of physical processes, 
for the longwave radiative cooling, was rewritten and modularized 
code will be considerably easier to maintain and modify. 


including that 
• The new 
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Response of Winter Forecasts made with the GLA 4th Order GCM 
to Changes in the Horizontal Grid Resolution 

J. Pfaendtner 


The GLA Fourth Order General Circulation Model was modified in 1984 to allow 
it to be run at a finer horizontal grid resolution. Previously the standard model 
versions integrated the primitive equations on a 4 degree latitude by 5 degree 
longitude grid. In this paper forecasts made at this resolution are compared 
to others which were made using the finer 2 degree latitude by 2 1/2 degree 
longitude grid resolution. 

Three winter cases are included with initial conditions from 15 December 
1978, 5 January 1979 and 21 January 1979. The 21 January 1979 initial conditions 
were, taken from a GLAS analysis (experiment 2728) which incorporated the GLAS 
temperature retrievals. The other two cases used European Center Analysis for 
initial conditions. The forecasts were all verified against European Center 
Analysis. 

Table 1 presents a summary for the three winter cases over North American 
and Europe. In general, the forecasts done on the finer grid were better out 
to day six. The first four days show significant reductions in the RMS errors. 
After day six, the 2 by 2 1/2 degree forecasts are worse than those done at 4 
by 5 degrees. A possible cause for this degredation is the stronger smoothing 
(8th order Shapiro filter rather than 16th order) used in the finer resolution 
model. The RMS errors for the individual forecasts are shown in tables 2 
through 5. A summer case (European Center initial condition from 15 June 1979) 
is also included in these tables. 

Table Is Percent Change in RMS Forecast Errors (2x2 1/2) vs. (4 x 5) 

Average of Three Winter Cases 

verification day 


1 2 3 4 5 6 7 8 9 10 


North Amer. 
pres . 

-6 

-16 

-9 

-6 

+2 

CM 

1 

+8 

-14 

+13 

-1 

Europe 

pres. 

-3 

-23 

-38 

-29 

+4 

-13 

+ 1 

+36 

+20 

+28 

North Amer. 

500 mb height 

-13 

-21 

-26 

-10 

-10 

-18 

+5 

-15 

+15 

+38 

Europe 

500 mb height 

-15 

-29 

-27 

-25 

-20 

-4 

+3 

+1 

-5 

-4 
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Table 2: RMS Sea Level Pressure Errors for North America (mb) 


verification day 

1 2 3 4 5 6 7 8 9 10 


1 

3553 15 Dec 78 


5.5 

5.8 

10.1 

13.7 



11.5 

12.1 

15.9 

3509 05 Jan 79 

m 

4.1 

5.8 

HI 

8.0 

8.1 




14.8 

3533 21 Jan 79 

4.1 

7.1 

H 

10.1 

m 

H 

8.0 

9.1 

m 

10.7 

Ave. of 3 winter cases 

H 

m 

7.0 

8.8 

9.5 

10.1 

9.6 

11.6 

12.0 

13.8 

3577 15 June 79 

2.8 

3.1 

HI 


5.0 

5.8 

8.1 

m 

4.1 


2x21/2 

cases 

3643 15 Dec 78 


H 

5.5 

10.4 






13.7 

3511 05 Jan 79 

2.7 

3.5 

H 

H 



n 



Hi 

3531 21 Jan 79 

3.5 

6.0 

8.1 

8.2 


6.8 



HI 


Ave. of 3 winter cases 

3.5 

m 

m 

8.3 

m 





13.7 

3727 15 June 79 

HI 

3.8 

4.8 

m 

7.0 

9.0 

8.4 

9.1 

9.3 

12.8 


Table 3: RMS Sea Level Pressure Errors for Europe (mb) 


verification day 

1 2 3 4 5 6 7 8 9 10 


4x5 

cases 

3553 15 Dec 78 


H 

m 

10.1 

11.7 

15.7 

18.3 

17.1 


16.6 

3509 05 Jan 79 

2.5 

5.2 

12.7 

11.5 


HI 



15.1 

17.4 

3533 21 Jan 79 

3.0 

4.8 

6.5 

10.0 

9.3 





15.7 

Ave. of 3 winter cases 

2.9 

4.8 

8.9 

10.5 

10.1 

11.4 

14.6 



16.6 

3577 15 June 79 

2.5 

4.2 

5.0 

5.2 

4.6 

5.9 

9.0 

9.1 

9.0 

8.9 

2x21/2 

cases 

3643 15 Dec 78 

3.4 

3.5 

5.3 

5.8 

10.8 

12.5 

17.4 

20.5 

21.4 

21.5 

3511 05 Jan 79 

2.3 

2.8 

5.5 

8.5 

12.1 

9.0 

11.9 

17.1 

16.2 



3531 21 Jan 79 

2.7 

4.8 

5.8 

7.9 i 

8.6 

8.2 



_ 

_ 

_ 

Ave. of 3 winter cases 

2.8 

3.7 

5.5 

7.4 ' 

10.5 

9.9 

14.7 

18.8 

18.8 

21.5 

3727 15 June 79 

2.4 

4.0 

4.2 

5.5 

4.7 

5.8 

7.9 

8.5 

11.9 

11.8 
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Table 4: RMS 500 mb Geopotential Height Errors tor North America (m) 

verification day 

1 2 3 4 5 6 7 8 9 10 


4x5 

cases 

3553 15 Dec 78 

41.7 

68.7 

91.0 

124.0 

154.0 

158.0 

126.0 

113.0 

125.0 

135.0 

3509 05 Jan 79 

24.2 

50.0 

60.0 

66.7 

82.6 

81.9 

98.9 

146.0 

145.0 

132.0 

3533 21 Jan 79 

41.1 

78.0 

110.0 

134.0 

134.0 

130.0 

107.0 

125.0 

124.0 

121.0 

Ave. of 3 winter cases 

35.7 

65.6 

87.1 

108.0 

124.0 

123.0 

111.0 

128.0 

131.0 

129.0 

3577 15 June 79 

25.2 

32.4 

46.8 

52.7 

54.8 

50.6 

62.4 

58.2 

50.9 

52.5 

2x21/2 

cases 

3643 15 Dec 78 

36.0 

54.3 

71.8 

125.0 

180.0 

168.0 

155.0 

122.0 

152.0 

178.0 

3511 05 Jan 79 

22.6 

38.2 

44.1 

62.7 

62.6 

60.9 

76.8 

96.4 

150.0 


3531 21 Jan 79 

34.2 

63.7 

1 

78.3 

104.0 

91.4 

74.3 





Ave. of 3 winter cases 

30.9 

52. 1 

64.7 

97.2 

111.0 

101.0 

116.0 

109.0 

151.0 

178.0 

3727 15 June 79 

21.5 

1 

27.5 

48.3 

58.5 

l 

66.9 

82.0 

82.1 

61.9 

64.4 

81.4 


Table 5: RMS 500 mb Geopotential Height Errors for Europe (m) 

verification day 

12 3 4 5 6 7 8 9 10 


4 x 5 
cases 

3553 15 Dec 78 



93.6 

124.0 

144.0 




in 


3509 05 Jan 79 



133.0 

146.0 

119.0 




ms 

RMra 




65.4 

91.8 

120.0 

111.0 

161.0 

182.0 

179.0 

196.0 

Ave. of 3 winter cases 

27.8 

52.5 

97.3 

121.0 

128.0 






3577 15 June 79 

23.4 

43.3 

65.3 

63.6 

63.5 





101.0 

2x21/2 

cases 

3643 15 Dec 78 

20.7 

41.8 

79.5 

92.4 

107.0 

143.0 

188.0 

211.0 

217.0 

201.0 

3511 05 Jan 79 

20.1 

31.8 

71.4 

89.7 

115.0 

105.0 

128.0 

174.0 

171.0 


3531 21 Jan 79 

30.4 

38.0 

61.3 

1 

88.5 

86.2 

85.2 





Ave. of 3 winter cases 

23.7 

37.2 

70.7 

90.2 

103.0 

111.0 

158.0 

193.0 

194.0 

201.0 

3727 15 June 79 

P 

23.9 

39.8 

61.8 

64.5 

64.6 

57.4 

80.0 

103.0 

151.0 

168.0 
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D. ATMOSPHERIC DYNAMICS AND DIAGNOSTIC 

STUDIES 




Comparison of Forecast and Observed Energetics 
W. E. Baker and Y. Brin 


The objective of this investigation is to perform diagnostic studies aimed 
at furthering the understanding of the atmospheric general circulation and 
providing useful insight into the NASA/Goddard analysis/forecast system as it 
continues to evolve* 

An energetics analysis scheme has been developed for this purpose to compare 
the observed kinetic energy balance over North America with that derived from 
forecast fields of the GLAS fourth order model for the 13-15 January 1979 
cyclone case* The major findings of that investigation are: 

1) The observed and predicted kinetic energy and eddy conversion are in good 
qualitative agreement, although the model eddy conversion tends to be 2 to 3 
times stronger than the observed values (see Figs. 1 and 2). The eddy conversion 
being stronger in the 12 h forecast than in observations may be due to a number 
of factors (e*g. an imbalance from the initial data, overestimation of the 
release of available potential energy from the physical parameterizations, 

etc*) which we plan to investigate* 

2) In agreement with previous studies of cyclonic disturbances (see Fig* 3) 
vertical profiles of kinetic energy generation and dissipation exhibit lower 
and upper tropospheric maxima in both the forecast and observations. 

3) An interesting time lag is noted in the observational analysis with the 
maximum in the observed kinetic energy occurring at 0000 GMT 14 January over 
the same region as the maximum eddy conversion 12 h earlier (compare Fig* la 
with Fig. 2a). 

Future work includes the examination of the forecast error over North 
America in terras of the energetics utilizing the limited-area, analysis scheme 
previously developed* The observed energetics will be compared with those from 
model forecasts in which the resolution, physics, or initial data have been modified. 


Reference 


Baker, W. E. , and Y* Brin, 1985: A comparison of observed and forecast energetics 
over North America. Quart. J* R* Meteor, Soc* , in press. 
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Fig, 1. Vertically integrated (surface to 100 mb) kinetic energy in 10 6 JnT 2 for 0000 GMT 14 January 1979. 

a. Observations (station analysis). 

b. 24h forecast. 
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Fig. 3. Vertical profiles of area-averaged kinetic energy 
generation and dissipation for 1200 GMT 13 January 
1979. The station values are denoted by a solid 
line and the 12h forecast by a dashed line. 
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Error Growth in Operational 10-day Forecasts 

E. Kalnay and A. Dalcher* 

Laboratory for Atmospheres 
NASA/Goddard Space Flight Center 
Greenbelt, Maryland, U.S.A. 


Lorenz (1981) estimated the rate of growth of forecast errors using as data 
base 100 consecutive 10-day operational forecasts from the European Centre for 
Medium Range Forecasting (ECMWF). He parameterized the growth of the difference 
between two model forecasts verifying the same day. We extended Lorenz* 
parameterization of error growth by including the effect of model deficiencies 
as well as their scale dependence . In this way we are able to fit real forecast 
error growth and determine the scale dependence of the derived parameters (rate 
of error growth, external error source due to model errors, and saturation 
levels). This work is a necessary step in the application of the Lagged Average 
Forecasting (LAF) technique (Hoffman and Kalnay, 1983, 1984) to the ECMWF 10-day 
forecasts (Dalcher et_al . , 1985). Like Boer (1984) and Savijarvi (1984), we 
found that for both the analyses and the forecasts, there is an approximate 
equipartition of energy between all the zonal wavenumbers m corresponding to 
the same total wavenumber n. 

In this work, we parameterize the growth of the error variance V by 
including both an external source of error and the saturation effects: 


— = a (V + S) (V„ - V) 
dt 


This equation has the solution 


( 1 ) 


v(t) = V.-H- - S _i_ 

1+p 1+p 


( 2 ) 


where p = ce a (V,,, +S)t, and the constant c = y(°) + S re i a t;ed to the initial 

Voo - V(o) 

error variance V(o). 


The solution (2) can also be expressed as 

V(t) = 1/2 (Voo + S)(l + tanh (In c + a (Voo + S)t)/2) - S 


(3) 


This parameterization work works remarkably well when all scales are 
combined together. Furthermore, the shape of the individual error growths for 
each total wavenumber suggests that the same parameterization can be applied to 
each total wavenumber. In Fig. 1 we apply the same function to three individual 
wavenumbers n, also showing an excellent fit, and indicating that the saturation 
error is attained at different times for different scales. 


*M/A COM Sigma Data Corp. at NASA/Laboratory for Atmospheres 
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Equation 1 is scaled for each wavenumber by dividing the error variance by 
the maximum observed value, which generally occurs at day 10. Some resulting 
parameters are presented in Fig. 2. In Fig. 2a we see that the scaled error 
growth rate a increases rapidly with wavenumber and is larger for the external 
errors (comparing forecasts with analyses) than for the internal errors 
(comparing two forecasts). 

The last result, Fig. 2b corresponds to the predictability time as a 
function of wavenumber, defined as the time needed to reach 95% of the error 
saturation value. The external error curve shows that for low wavenumbers 
(n K 9), there is still predictability in the ECMWF model after 10 days, whereas 

larger wavenumbers become saturated at increasingly shorter times. The curve 
corresponding to internal errors suggests that with a perfect model, predicta- 
bility would be significantly larger, reaching over 3 weeks for n < 5. 
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? ig. 1. Fit of the error variance growth for waves with total wave number 
n =* 3, 7 and 15. The data values are marked by dots. The curve 
is the fitted solution and is extrapolated up to 30 days. 
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Mechanistic Experiments to Determine the Origin of Southern Hemisphere 

Stationary Waves 

Eugenia Kalnay and Kingtse C. Mo 
Laboratory for Atmospheres 
NASA/Goddard Space Flight Center 
Greenbelt, MD 20771 


Kalnay and Paegle (1983) reported the existence and structure of large 
amplitude stationary Rossby waves with zonal wavenumber ~ 7 between 20°S and 
40 S over and in the lee of South America* (Fig.l) These waves, present during 
January 1979, disappeared in February 1979. Kalnay and Paegle (1983) discussed 
their possible origin, and concluded that the waves were not orographically 
forced because they corresponded to a ridge in the lee of the Andes, and were 
not forced by observed sea surface temperature anomalies because of their 
space and time phase relationship. They suggested that the waves were probably 
forced by tropical heating. 

In this paper we perform two 15 day forecast experiments with the GLAS 
^ou^th Order General Circulation Model, and initial conditions corresponding to 
5 January and 4 February 1979. These forecasts reproduce reasonably well the 
presence of the January wave. (Fig. 2) and their absence in February. Several 
mechanistic experiments to determine the origin of the waves are then performed. 

A "No Andes” forecast shows that the waves do indeed exist independently of oro- 
graphic forcing. (Fig. 3) Several experiments modifying the coefficient of 
latent heat lead to the conclusion that tropical heating is important in the 
maintenance of the waves (Fig. 4). Furthermore, the convection in the subtropical 
waves themselves is important in sustaining their amplitude and phase, and the 
Walker type of circulation associated with the SPCZ is also a contributor to 
the maintenance of the South American waves. These results confirm the existence 
of a relationship between the occurrence of a strong South Pacific Convergence 
Zone (SPCZ), (Fig. 5) somewhat eastward from its climatological position, and 
the strong "South Atlantic Convergence Zone" (Fig. 6) observed in outgoing long 
wave radiation maps* 

Reference 


Kalnay, E., and J. Paegle 1983: Large amplitude stationary wind in the Southern 
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?ig.l 15 day averaged meridional velocity Fig. 2 Same as Fig.l but for the January 

at 200 mb from the GLAS analysis, control forecast, 

from 5 to Jan 20, contour interval 
5 m/s. 



Fig. 3 Same as Fig. 2 but for the "no Andes" Fig. 4 Same as Fig. 2 but for the reduced 

experiment. tropical heating experiment. 


tort rcm * | mm mm* 



Fig. 5 Same as Fig. 2 but for the "Suppressed Fig. 6 Same as Fig. 2 but for the "Suppressed 
Pacific heating experiment. Atlantic heating" experiment. 
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Energy Sources of the Dominant Frequency Dependent 
3-Dimensional Atmospheric Modes 

Siegfried Schubert 


As part of an on-going study of atmospheric variability we present here the 
energy sources and sinks associated with the zonally asymmetric winter mean flow. 
The observational work of Blackmon et al . , (1984) showed distinctly different 
horizontal structures for the long, intermediate and short time scale atmospheric 
variations. In this study we look into the 3-diraensional structure of the 
fluctuations and examine the relative roles of barotropic and baroclinic terms. 

In the context of a two-level quasi-geostrophic model the energy tendency 
terms are 

• • • 2 s 2 s 

<K X + + P> = 1/4-rr j {V i|> (kx Vi(i • Vip ) + V x (k x ViJ> • Vt ) 

2 s 2 s 2 s 

+ V (k x Vt • Vt ) + V t (k x Vt • Vi|> ) + r x (k x ViJ> *Vt) } d< (1) 

where K T and are the vertical mean flow and shear flow kinetic energy (KE) 
associated with the anomalies, respectively and the superscript (s) denotes the 
winter mean flow. The available potential energy (APE) associated with the 
anomalies is P = r^x^/2 where r is the non-dimensional inverse Rossby deformation 
radius. The anomalies are defined as deviations from the winter mean flow where 
x and V are one-half the vertical sum and difference of the 200 mb and 700 mb 
stream function fields, respectively. The data consists of 10 years of Northern 
Hemisphere winter stream function fields where, in practice, 'F is approximately 
by the 500 mb field. The first four terms on the RHS of (1) combine to give the 
total barotropic conversion, while the fifth terra represents the conversion of 
the time mean flow APE to the anomaly APE. 

In order to isolate the relative roles of these terms for the wide range of 
frequencies which exist in the atmosphere, five different filters were applied to 
the anomalies. Filters A— E were designed to retain fluctuations longer than 45 
days, 20—45 days, 10—20 days, 6—10 days, and 2.5—6 days, respectively. For each 
of the five frequency bands the data was expanded in empirical orthogonal functions 
and only the first five EOFs were retained. 

Figure 1 compares the time averaged total KE conversion term and the APE 
conversion term for the first five EOFs in each frequency band for a deformation 
radius typical of middle latitudes (1000 km). Positive values indicate a source 
(sink) for the anomalies (mean flow). These results show a dramatic shift in the 
nature of the source terms from a predominantly barotropic source for the lowest 
frequency modes to a predominantly baroclinic source for the synoptic time scales. 
Further details of these conversion terras including their spatial distributions 
may be found in Schubert, 1985. 
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Principal Component 

Figure 1. The time-average energy sources associated with the winter mean flow 

for the first five PCs for filters A through E. PCs 1-5 correspond to filter 
A (lowest frequencies), PCs 6-10 correspond to filter B, etc. See text for 
filter characteristics. The ordinate is the energy tendency divided by the 
variance of the associated PC. 
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A Comparison of the Bounded Derivative and 
the Normal Mode Initialization Methods Using Real Data 

F. H. M. Semazzi and I. M. Navon 


Browning et al . (1980) proposed an initialization method called the bounded 
derivative method (BDI). They used analytical data to test the new method. 

Kasahara (1982) theoretically demonstrated the equivalence between BDI and the 
well known nonlinear normal mode intialization method (NMI). The purposes of 
this study are the extension of the application of BDI to real data and comparison 
with NMI. The unbalanced initial state (UBD) is data of January 1979 00Z which 
were interpolated from the adjacent sigma levels of the GLAS GCM to the 300 mb 
surface. We use the global barotropic model described by Takacs and Balgovind 
(1983). Orographic forcing is explicitly included in the model. 

Many comparisons are performed between various quantities. However, we only 
present a comparison of the time evolution at two grid points A(50° S, 90° E) and 
B(10° S, 20° E) which represent low and middle latitude locations. To facilitate 
a more complete comparison an initialization experiment based on the classical 
balance equation (CBE) was also included. Figs, la, lb and lc show the time 
evolution of the height corresponding to UBD, CBE, BDI, and NMI at coordinate A. 

The forecast starting from UBD suffers from contamination by high frequency 
oscillations with a dominant period of about 6 hours. Regarding CBE, BDI and NMI 
the time evolution is smooth. The classical balance equation performs well since 
the influence of orography is small at A. The plots for the time evolution for 
coordinate B are presented in Figs. 2a, 2b and 2c. Notice the balance equation 
performs relatively worse than the NMI and BDI. This is mainly because the 
balance equation is not capable of handling the gravity waves produced by orography 
In this case these waves are excited by the high elevation over equatorial Africa. 
It is also instructive to note that at low latitudes the amplitude of gravity 
waves can be larger than that of synoptic motions. However, over higher latitudes 
the relative importance of gravitational motions is smaller. This indicates the 
importance of performing suitable initialization over the low latitude regions 
of the earth, particularly for purposes of short range numerical weather predic- 
tion. 
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The Vertical Structure of Global Rotational Normal Modes 

David M. Straus 


The relevance of global rotational normal modes to the behavior of the 
atmosphere has been a concern of meteorologists for nearly half a century. 
Although the mathematical properties of these fundamental linear solutions 
have been explored at various levels of approximations over the years, it is 
only recently that the availability of reliable global data has made meaningful 
global observational studies possible. 

In one such recent study, Lindzen et al . (1984) examined the amplitude 
and phase evolution of Hough mode projections at 500 mb. The Hough modes 
represent the simplest normal mode approximation available, and correspond to 
the neutral eigenfunctions of a shallow water fluid with no mean zonal flow. 

In this study, it was found that when the observed amplitude of a rotational 
Hough mode was large, it tended to propagate at the phase speed of a normal 
mode in the presence of mean 500 mb winds, giving a strong indication that 
normal modes are of relevance to the atmosphere. 

Continuing along these lines, we have explored the vertical structure of 
the approximately defined rotational normal modes by projecting observed data 
(from the ECMWF's FGGE analyses) onto Hough functions at levels other than 
500 mb. As in Lindzen et al . (1984), the stationary and eastward propagating 
components were filtered out at each level. The evolution of the amplitude in 
time throughout the troposphere is shown for several modes during summer and 
winter in Figures 1-4. The episodic nature of the normal modes at 500 mb found 
by Lindzen et al . is evident at other levels, and it continues to be true that 
large amplitude is associated with phase propagation at the theoretical phase 
speed (corrected for the presence of mean zonal winds). 

What is suprising about the vertical structure depicted in these figures 
is that bhe level of maximum amplitude is 300 mb, so that the amplitude decreases 
with altitude in the upper troposphere. This is contrary to the predictions 
of both the simple theory and more complete theoretical-numerical determinations 
of normal mode vertical structure that the amplitude increases throughout the 
troposphere and even, in some cases, the stratosphere. 

This anomalous distribution of amplitude is potentially of great interest, 
for it may indicate either that some vital ingredient has been left out of 
standard normal mode calculations, or that the observed modes’ phase progression 
at the normal mode speed is fortuitous and the identification as normal modes 
incorrect. (These two possibilities are not necessarily distinct.) A third, and 
perhaps more likely alternative, is that the vertical structure observed in the 
figures is the result of using Hough functions, rather than eigenfunctions 
which take into account the presence of the mean zonal wind, which itself tends 
to be a maximum above 300 mb. 

Current research centers on projecting data from each level onto shallow 
water eigenfunctions which take into account the mean zonal wind at that level in 
order to clarify the question of vertical structure. This is of importance for 
future work on the role of these modes in numerical forecasts. 
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WINTER 
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Fig. 1 Evolution of amplitude for the 
(2, 4) (wavenumber 2, first 
antisymmetric) mode as a funct- 
ion of level for the period of 
Dec. 31, 1978 to Jan. 28, 1979. 
The amplitude represents the 
value of the geopotential 
height in meters. The Hough 
mode projections were carried 
out at 700, 500, 300 and 100 mb. 
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Fig. 2 Evolution of amplitude for the 
(2, 5) (wavenumber 2, second 
symmetric) mode for the period 
of Dec. 31, 1978 to Jan. 28, 
1979). 
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Fig. 3 Evolution of amplitude for the 
(2, 5) mode for the period of 
June 1 to June 30, 1979. 


SUMMER 
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Fig. 4 Evolution of amplitude for the 
(3, 5) (wavenumber 3, first 
antisymmetric) mode for the 
period of June 30 to July 30, 
1979. 
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II. CLIMATE/OCEAN-AIR INTERACTIONS 


A. DATA ANALYSIS 




Quasi stationary states in the Southern Hemisphere 


Kingtse C. Mo 


Pattern correlations between daily anomalies (Horel 1984) have been used 
to study the persistency of the Southern Hemisphere circulations. The data set 
consists of the daily Australian analyses of 500 mb heights and sea level 
pressure for the period from 1972 and 1983. Compared to the Northern Hemisphere, 
the pattern correlations are much smaller and more variable in the Southern 
Hemisphere. The mean one day lag is only 0.57, (0.82 for the NH). The corre- 
lations increase significantly for the filtered anomalies which consist of the 
long planetary waves from 0 to 4. 

Subjective criteria based on the pattern correlations are used to select 
the quasi stationary events. A series of 5 or more daily maps is quasi stationary 
if the pattern correlations between all consecutive pairs of maps in this time 
series are larger or equal to 0.5. 

In winter, events can be classified in terms of wave numbers. Wave 3 and 
4 are by far the most dominant waves. There are total 24 events in winter. 

There are 14 events which have large 3 amplitude. 8 out of these 14 events 
have fixed phase locations. The composite (Fig. 1) of the total anomalies for 
these 8 events shows lows at 50S, 115E, 60S, 0W and 55S, 140W. The same wave 3 
signal has been found in the winter teleconnection pattern using monthly mean 
data. (Mo and White, 1984). 

There are 3 events which also have a large wave 3 amplitude but opposite 
phase locations. The composite of the total anomalies of these 3 events is 
given in Fig. 2. The reversal of positions of highs and lows may suggest the 
superresonant and subresonant states described by the multiple equilibria 
theory. (Charney and Straus 1980). 

In summer, these are only 12 events. All of them fall in months where the 
Southern Oscillation reaches extremes. 
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Trends in the Southern Hemisphere 
Kingtse Mo 


and 

Harry van Loon 


I. Introduction 

The station data of monthly mean sea level pressure, surface air tempera- 
ture as well as monthly mean upper air rawinsonde observations were used to 
study the trends in the southern winter (June, July and August) and summer 
(December, January and February). This study covers the period 1951 to 1981. 

Most of the surface data have continuous 31 year records, however, for some 
high latitude stations, both surface and upper air records started in 1957 or 
even later. 

The normalized anomaly for a given variable Z is defined as the departure 
from the 1951-1981 grand mean divided by the standard deviation during the same 
period. To eliminate the large fluctuations in the data, the 1-3-1 smoothing 
has been applied to the time series of normalized anomalies. 

II. The Sea level pressure 

Figure 1 plots the smoothed normalized anomalies of sea level pressure for 
31 stations during the winter season. The thick line gives the parabola fit 
for unsmoothed anomalies. The trends are highly regionally dependent. Except 
for stations at high latitudes, the mean sea level pressure fell at stations 
on the Indian and Atlantic side of the hemisphere from the 1950 f s to the 1960 T s 
but it rose over the other half of the hemisphere. These trends in general 
reversed from the 1960 f s to the 1970 f s. 

At higher latitudes, e.g. Campbell and Macquarie on the Pacific side, 
sea level pressure increased during these 30 years. On the opposite side, 
e.g. Ushuaia and Stanley, the pressure fell. At Halley Bay and S.A.N.A.E. the 
sea level pressure was higher during the 1960 ! s than the 1970 f s. 

These trends affected the amplitude of the half yearly wave in the midlat- 
itudes (Mo and van Loon, 1984, Swanson and Trenberth, 1981). They also influenced 
on amplitude of the planetary waves. 

III. Other Points 

In this study, we document the trends of monthly mean temperatures and 
heights from the surface to 200 mb level. Other conclusions are: 1. Trends 
of all these elements are highly regionally dependent. The mean temperature 
of the earth and its changes can not be adequately described if the selected 

stations are not evenly distributed over the globe. 2. The trends of the 
height fields are equivalent barotropic. 3. There are no major differences 
in trend between stations on land or sea (Chen 1982). 4. A linear trend 

indicates a large increase in heights throughout Antarctica from 1957 to 1981 


93 



at all levels. This shows an increase of the total mass in the atmosphere 
over that region (Kraus 1977). 5. There was a notable, warming over Antarctica 

from the 1960's to the 1970's except in a small area from 90W to 120W at all 
upper levels. (Hanson et al .. 1981). 
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Fig. 1 The smoothed normalized anomalies (thin line) of sea level pressure in 
winter, the thick line gives the parabala fit for unsmoothed anomalies. 




The Seasonal Cycle of Storminess as Measured by Band-Pass Fluctuations 

David M. Straus 


An important prerequisite to the successful modeling of climate variations 
is the ability to simulate the observed seasonal cycle in both the time and 
space averaged fields and in the fluctuations of these fields. With the advent 
of high speed vector computers, a number of multi-year simulations with differen- 
GCMs has been carried out with the seasonal cycle in boundary conditions (solar 
insolation/ sea surface temperature/ etc.) imposed. The cottparison of the 
seasonal cycles produced by the GCMs with the observations should be a priority 
item on the agenda of modeling groups. 

This endeavor will require the compiling of appropriate observational 
statistics. The seasonal variation of the mean fields is well documented and 
several climatologies of local time variances during winter and summer are 
available. However, good statistics for the transition seasons of spring and 
fall are harder to obtain. We are undertaking the compiling of climatologies 
of local time variances and covariances for all four seasons from a 15 year 
dataset consisting of NMC operational (Northern Hemisphere) analyses from May 
1963 through December 1977. The emphasis is to be on the clear depiction of 
the seasonal change of the ( co) variances. 

In this note we present a sample of such statistics, namely the seasonal 
cycle of baroclinic storms, as represented by band-pass filtered geopotential 
height variances at 850 mb. The particular filter used is that suggested by 
Blackmon and White (1982), and retains periods of approximately 2.5 to 10 
days. The time series of height (at each grid point) were filtered by removing 
the annual and semiannual cycles for that point, and by removing zonal wavenumber 
higher than 20. The band-pass filter was then applied. The height variances 
of the filtered fields were then computed for each winter season (December- 
February), each spring season (March-May), each summer season (June-August) 
and each fall season ( September-November ) . These variances were then averaged 
by season. Figures 1(a) -(d) show maps of the standard deviation. These figures 
show clearly the seasonal cycle of bandpass fluctuations. The strong maxima 
over the Atlantic and Pacific Oceans are present in all seasons, with a signifies 
weakening (and slight northward shift) only in the summer. The major seasonal 
variation is thus seen to consist mostly of a summertime weakening and shift? 
spring and fall appear nearly identical to winter. The corresponding results 
at 500 mb (not shown) are similar, but with the stormtrack variance being 
slightly larger in spring and fall compared to winter. 
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The Contrasting Northern Hemisphere Winters of 1980-81 and 1981-82 

Glenn H. White 


Seasonally averaged statistics calculated from daily operational analyses 
by the European Centre for Medium Range Weather Forecasts (ECMWF) for December 
1980-February 1981 (winter 1) and December 1981-February 1982 (winter 2) are 
used to contrast the two seasonal mean circulation patterns present in the 
Northern Hemisphere and to investigate possible causes for the differences. 

The differences between winter 1 and winter 2 resemble patterns associated 
with the Southern Oscillation (White, 1984) and may be forced by a weak warm 
El Nino-Southern Oscillation (ENSO) event in 1980. 

Low-level diabatic heating is much stronger over the North Atlantic in 
winter 2 (Fig. 1), where changes in sea surface temperatures imply a much 
greater loss of heat by the ocean to the atmosphere in winter 2. Changes in 
the flow pattern over the United States and the North Atlantic from winter 1 
to winter 2 resemble the response to warm sea surface temperature anomalies 
found by Phillips (1982) in a nonlinear model of a baroclinically unstable 
atmosphere. A large change of opposite sign in diabatic heating, strongest at 
500 mb, is observed over Japan. 

Zonally averaged and three-dimensional Eliassen-Palm fluxes (Plumb, 1984) 
for the two winters (Fig. 2) suggest that stationary waves are mainly forced 
near the surface in midlatitudes. The forcing appears to be significantly less 
in winter 2, particularly near the south coast of Alaska where the low-level 
seasonal mean flow over orography was much weaker in winter 2. Since changes 
in the low-level flow over the Pacific appear related to changes in equatorial 
sea surface temperatures and precipitation in the Pacific (Fritz (1984)), a 
weak ENSO event in 1980 may have produced a delayed response in the circulation 
over the Pacific in winter 1 that resulted in unusually strong orographic 
forcing over Alaska and western Canada. 
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Fig. 1 Change in diabatic heating in the 1000-775 mb layer from 

winter 1 to winter 2. Contour interval 15 watts/m^, dashed 
contours represent cooling. 
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winter 1 (Fig. 2a) and winter 2 (Fig. 2b) 




On the Global Distribution of Three-dimensional 
Eliassen-Palm Fluxes by Stationary Waves 

Glenn H. White 

New insights into the propagation of atmospheric waves have recently been 
gained from the use of zonally averaged Eliassen-Palm (EP) fluxes (Edmon et al. , 
1980). If wave dynamics can be described by linear planetary-wave theory and if 
the poleward gradient of quasi-geostrophic potential vorticity is positive, EP 
fluxes measure the net propagation of waves and can be used to infer the forcing 
of waves. A three-dimensional EP flux for stationary waves has recently been 
derived by Plumb (1984): 

F = (_P) cos <J> f(v*2 - _L <Kv '? 1 ^ + f-u*v* + _J_ 3(u T ) $ 

~ Po 2f 3x 2f 3x J 


+ I fv*T* - 3(?*r* ) ) 

S 2f 3x 

where ( ) represents a time-mean, ( )* the departure from the zonal mean, u and 
v the zonal and meridional winds, T temperature, $ geopotential, <j> latitude, f 
the Coriolis parameter, S static stability (allowed to vary only with pressure 
p), and p Q = 1000 mb. Plumb (1984) presented observational results for the 
extratropical Northern Hemisphere wintertime climatological circulation based 
on statistics calculated from NMC operational analyses for 1965-1976. 

This study investigated the global distribution of three-dimensional EP 
fluxes (Fig. 1) by using statistics for Dec. 1980-Feb. 1981 (Figs, la, b) and 
June-August 1981 (Figs. lc,d) calculated from operational analyses by the 
European Centre for Medium-range Weather Forecasts. Results at 1000 (Figs. 
la,c) and 150 mb (Figs. lb,d) are shown. During the Northern Hemisphere (NH) 
winter strong upward fluxes can be seen at 1000 mb (Fig. la) in the lee of the 
Himalayas and Rockies, windward of the Canadian Rockies in the east Pacific, 
to the north of the Atlantic stormtrack and over the west Pacific in a region 
of strong land-sea thermal contrast. Upward fluxes also appear near the west 
coasts of the Southern Hemisphere (SH) subtropical continents regions of strong 
land-sea thermal contrast. At 150 mb (Fig. lb) the strongest upward flux now 
occurs south of Alaska. The horizontal fluxes imply poleward propagation from 
the convection over Indonesia to 20°N, a region where the assumptions underlying 
EP fluxes are not well-met. Equatorward propagation dominates the NH midlatitudes 
except over the Pacific. 

During June-August 1981 at 1000 mb (Fig. lc) upward fluxes appear off the 
western coasts of the subtropical continents in both hemispheres but are parti- 
cularly strong in the NH where strong land-sea thermal contrast is found. At 
150 mb (Fig. Id) strong fluxes from the Indian monsoon into the SH can be seen 
where easterly winds prevail, which theory and numerical models predict would 
prevent propagation (Webster and Holton, 1982). Strong fluxes also occur at 
30°N over the eastern Mediterranean, a region of strong sinking, and over the 
Pacific and the southeastern U.S., regions of rising motion (White, 1983). 

The fluxes shown in Fig. 1 offer tantalizing hints of the sources of 
stationary waves, but also show several puzzling features and a rather cavalier 
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disregard of regions of easterly wind. The physical meaning and interpretation 
of three-dimensional EP-fluxes is not yet clear. 
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1. 3D Eliassen-Palm fluxes during Dec. 1980-Feb. 1981 at a) 1000 and b) 
150 mb and during June-August 1981 at c) 1000 and d) 150 mb. Solid 
contours denote upward fluxes and dashed contours downward fluxes; 
contour Interval at 1000 mb is ten times that at 150 mb. 











Transient eddies in the UCLA GCM 


Glenn H. White 


The simulation of transient eddies in the nine-level UCLA general circulation 
model (GCM) has been examined and compared to observations/ with emphasis on 
the Northern Hemisphere winter. Qualitatively, the UCLA GCM reproduces many 
features of the observed circulation and the relationship between the time-mean 
flow and transient eddies; however, the magnitudes of transient eddies in the 
UCLA GCM, particularly at frequencies lower than those associated with baroclinic 
instability, appear to be much less than those observed. 

Fig. la shows the wintertime mean 500 mb height field for one simulated 
winter. The Asian trough appears broader than in the observed climatology (Lau 
et_al . , 1981) while the ridge over western Canada is well to the east of its 
climatological position. The North American trough is much weaker than observed, 
while an unusually strong ridge occupies the east Atlantic.. Such anomalous 
features can be observed in seasonal mean fields for individual winters (Lau 
et al . , 1981); however, the other two winters in the simulation examined show 
very similar anomalous features, implying a departure of the climatology of the 
UCLA GCM from observations. 

The patterns of low-level transient eddy heat fluxes were similar to 
observations and acted to dissipate the low-level seasonal-mean temperature 
field, as observed (Lau et al . , 1981). The simulated shape of transient 
disturbances also appeared similar to observations, with eddies of baroclinic 
time scales displaying a north-south elongation and disturbances of longer time 
scales an east-west elongation (Hoskins et al . , 1983). Baroclinic eddies in 
the simulation dissipated when they encountered a longer-lived ridge, a pattern 
also observed (Hoskins et al . , 1983). 

The simulated transient eddy kinetic energy at 300 mb is displayed in 
Fig. lb for eddies of all time scales shorter than a season, in Fig. lc for 
time scales associated with baroclinic instability, and in Fig. Id for time 
scales longer than in Fig. lc. (The actual filters used were the "Lorenz" 
filters discussed on pg. 7 of Lau et al . , 1981.) The pattern in Fig. lb is 
somewhat similar to that observed at 250 mb by White (1983, pg. 10); however 
the maximum values in Fig. lb are only a third of the maximum values found by 
White (1983) and are located 10° further north. The maximum values for eddies 
of baroclinic timescales (Fig. 1c) are three-fourths of those observed (White, 
1983), suggesting, as Fig. Id shows, that the UCLA GCM is deficient in low- 
frequency variability, a shortcoming found in other GCM 1 s (Straus and Shukla, 

1981). The poleward shift of the simulated storm tracks relative to observations 
may reflect a similar shift in the region where the seasonal mean flow is 
baroclinically unstable, particularly in the Pacific. 
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Fig. 1 UCLA GCM simulations of wintertime mean a) 500 mb geopotential height 
and 300 mb transient eddy kinetic energy, b) unfiltered, c) high-pass 
filtered, and d) low-pass filtered. 





On Fofonoff's Mode 


Lee - Or Merkine 

Department of Mathematics Technion Israel Institute of Technology 
Kingtse Mo and Eugenia Kalnay 


The classical barotropic ocean model of Veronis (1966) is re-examined to 
verify the possible existence of Fofonoff's (1954) steady free inertial mode by 
considering the balance 


J ( n,RV 2 n + y)=0 (1) 

where J is the Jacobian, p is the streamf unction of the horizontal velocity, 
and R is a parameter of nonlinearity. A linear functional relationship between 
the potential vorticity and the streamf unction and inertial boundary layers 
along the ocean basin was suggested by Fofonoff. 

R V 2 n + y = G (p) 

where G (n) is an unknown function 

As in Veronis 1 work, we assume a square basin with the consequent boundary 
conditions. 


n = 0 on x = 0, tt and y = 0, tt 
and a wind stress given by 

*>• 

k • V x t = - sin x sin y 

We consider a very strong beta effect such that R = 10" 6 and e = 10“ 4 where 
s measures the strength of the dissipation. This set of parameters puts the system 
In the range where the Fofonoff’s hypothesis holds. 

However the numerical integration shows that the system is barotropically 
instable. The typical instantaneous field is shown in Fig. 1. The statistical 
equilibrium state n is depicted in Fig. 2. Although it resembles Fofonoff’s 
solution, it is maintained by the Reynolds stress field induced by the eddies in 
idditi on to the wind stress and dissipation. The plot of n versus R V 2 n +y 
Ls given in Fig. 3. It shows that the Fofonoff’s linear relation between the 
potential vorticity and the streamfunction does not hold uniformly over the 
inti re basin. The above results show that the Fofonoff’s solution can not be 
realized physically. 
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Analysis of a Simulated Cloud Climatology 
D. Randall, Harshvardhan, D. Short, and T. Corsetti 


Cloudiness simulations with the UCLA/Goddard GCM have been analyzed by 
comparison of simulated cloud and radiation statistics with corresponding 
observations. A description of the model formulation is given by Suarez et al . 
(1983), and some results are discussed by Randall et al . (1985). 

The global mean cloudiness is near 50%, in agreement with observations. 
However, comparison with the satellite observations of Susskind et al . (1983) 
shows that the model produces too many low-level clouds at high latitudes and 
too few high clouds in the tropics. The Intertropical Convergence Zone is 
hardly visible in the simulated cloudiness maps, although it is quite apparent 
in the simulated precipitation field. In the subtropics, the marine subtropical 
stratocumulus regimes are correctly positioned, but the cloud amounts are 
seriously underpredicted. The Arctic summer stratus is not well simulated. 

Over the Northern Hemisphere continents in winter the model produces an 
overabundance of low clouds, although these are geometrically and optically very 
thin, and might well be missed in the observations. 

The simulated outgoing longwave radiation at the top of the atmosphere is 
too low in middle and high latitudes, and much too high in regions of deep 
convection (Fig. 1). The standard deviation of daily mean outgoing longwave 
radiation is overpredicted in regions of deep convection. These results suggest 
that cirrus clouds are too infrequent and too variable in the model results. 

Seasonal variations of outgoing longwave radiation in the tropics and 
subtropics, as measured by the phase of the first annual harmonic, are fairly 
realistic in regions of observed low-level cloudiness, and quite unrealistic in 
regions of observed high-altitude cloudiness. 

In summary, the most serious deficiency of the current model's cloudiness 
simulation is its gross underprediction of the cirrus cloudiness associated 
with deep convection. We are currently working to remedy this problem by 
incorporating a prognostic variable for cirrus cloud water. 
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TERRESTRIAL RADIATION AT THE TOP OF THE ATMOSPHERE (W m -2 ) 
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Figure 1: Model simulated (right) and NIMBUS 7 - observed (left) terrestrial radiation at the 

top of the atmosphere for January (top) and July (bottom)* The contour interval is 
? S V - 











The Moist Available Energy 
of a Conditionally Unstable Atmosphere 

D. Randall and T. Corsetti 


Lorenz (1978, 1979) defined the moist available energy (MAE) of a given 
state of the global atmosphere as the difference between the moist enthalpy of 
the given state and that of a hypothetical '‘reference" state for which the 
moist enthalpy is minimized by reversible adiabatic rearrangement of the 
air parcels. Within the reference state, the pressure and density are horizontally 
stratified, hydrostatic equilibrium holds, and the stratification is statically 
stable. However, a given horizontal surface in the reference state may contain 
both a saturated cloudy region and an unsaturated clear region. 

For an atmosphere containing no moisture, the MAE is identical to the 
"available potential energy" or "dry available energy" discussed by Lorenz 
(1955). In general, however, the MAE of a humid atmosphere is greater than 
its dry available energy. For a realistic atmosphere which is everywhere 
statically stable (or neutral) in the dry sense, but contains regions of condi- 
tional instability, the moist available energy arises in part from horizontal 
temperature gradients and in part from the conditional instability. 

Consider a hypothetical horizontally homogeneous atmosphere for which the 
temperature and moisture soundings are typical of conditionally unstable regions 
in the tropics. For such an atmosphere the MAE is entirely due to the presence 
of conditional instability, and is a measure of the degree of conditional 
instability. How is the MAE of this atmosphere related to the cloud work 
function of Arakawa and Schubert (1974), or to other conventional measures 
of conditional instability? The cloud work function is a measure of the energy 
available to a particular type of circulation, represented by a very simple 
model of a cumulus cloud. The simple cloud model is actually used in the 
computation of the cloud work function. The MAE, on the other hand, is the max- 
imum energy available to any circulation, and its definition does not involve a 
particular cloud model. It is therefore more general than the cloud work func- 
tion, and must exceed the cloud work function in most cases. 

Lorenz (1978, 1979) suggested graphical and numerical algorithms for 
computing the MAE for a given state of the atmosphere. Both algorithms encounter 
some difficulty when the given state is conditionally unstable. We have devised 
two simple methods for computing the MAE for cases in which the total number of 
parcels is fairly small. The first is a "swapping" routine that finds the 
reference state by exchanging parcels whenever the swap reduces the total 
enthalpy. However, we can prove by contrived example that the swapping algorithm 
does not always find the true reference state. A second "brute force" routine 
actually examines all possible arrangements of the parcels and chooses the one 
with the smallest moist enthalpy. At present, both routines allow only horizon- 
tally homogeneous reference states. 

In an exploratory study, we have used the swapping routine to find MAE and 
reference state for a time-sequence of GATE soundings, all of which are condi- 
tionally unstable. In spot checks, the brute-force routine always agreed with 
the swapper. Results show that in passing from the given state to the reference 
state the lowest three hundred mb of the atmosphere typically exchanges places 
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with the 300 mb above it. The reference state then contains a dry warm layer 
near the surface, and a supersaturated but cool layer in the middle troposphere. 
The upper portion of the sounding is not modified in passing from the given 
state to the reference state. 

We are continuing this work by comparing the MAE with the cloud work 
function, by investigating the rates at which MAE is generated by surface 
fluxes and large-scale vertical motion, and by investigating the possibility 
of horizontally inhomogeneous reference states. 
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Preliminary GCM Results with a New Radiation Parameterization 
D. Randall, Harshvardhan, R. Davies, and T. Corse tti 


We have developed new parameterizations of solar and terrestrial radiation 
and they have been tested UCLA/Goddard GCM, The solar radiation parameterization 
is based on that of Lacis and Hansen (1974), but with zenith-angle-dependent 
surface albedoes, and revised treatments of cloudiness* The terrestrial radia- 
tion parameterization is based on the work of Chou 1984) for water vapor, Chou 
and Ihng (1983) for carbon dioxide, and Rogers (1968) for ozone, with a new 
parameterization of the effects of clouds. Some details are given by Harshvardhan 
and Corsetti (1984). 

The codes have been designed for efficient execution on a vector processor, 
and we find that, for all gridpoints in a 4 x 5 degree 9-level GCM, one pass 
consumes about 2 cpu seconds on Goddard's two-pipe CDC Cyber 205. This excellent 
computational speed allows us to do a full radiation calculation once each 
simulated hour. 

We have compared results obtained with the new parameterizations to those 
obtained with the earlier parameterization described by Schlessinger (1976). 
Several dramatic improvements have come to light. For the most part these 
are related to the fact that the new terrestrial radiation pararaerization includes 
the effects of the water vapor continuum, while the earlier parameterization 
does not. In the moist tropical planetary boundary layer, continuum emission 
leads to much stronger cooling of the FBL over the oceans. Over land, however, 
the cooling of the FBL is significantly reduced. The latter, somewhat para- 
doxical result is due to the strong diurnal cycle of the continental FBL. At 
night the shallow continental PBL is overlain by a moist layer created by 
mixing during the previous afternoon. The moist upper layer acts as a radiative 
blanket, reducing the time-averaged radiative cooling of the continental PBL. 

The stronger cooling of the marine tropical PBL leads to an increase in the 
surface sensible heat flux, and an increase in the relative humidity at the PBL 
top. These are both improvements in the realism of the model results. 

We have corrected an error which caused cumulus anvils to be ignored, with 
radiation calculation. As shown in Fig. 1, this has resulted in a simulated 
diurnal cycle of precipitation over the oceans, in agreement with the observation 
of Gray and Jacobsen (1977). 

In summary, changes to the model's radiation parameterization have led to 
significant improvements in the simulated climate. 
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Difference from Time-Averaged Total Precipitation, mm day -1 



Local Time 


Figure 1. Simulated departure of the total precipitation rate from the time 
mean, as a function local time of day, composited for 24 ocean 
grid points. 
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Multiple Equilibria of the Barotropic Vorticity Equation on a Sphere 

M. J. Suarez 

In their landmark paper on multiple-equilibria, Charney and Devore (1979) 
(CD) argued that in the presence of topography a forced zonal flow can satisfy 
the zonal momentum balance in more than one way. For a strong zonal flow the 
stationary wave is small and results in little form drag on the zonal flow so 
that it can be at equilibrium near its forced value. A second state can occur 
provided the mountain is sufficiently high, the wave damping sufficiently 
small, and the forced wave is contained in a resonant cavity. In this state 
the wave amplitude is large, and the form drag on the zonal flow can balance a 
large departure from the forcing. Such a state is possible only if the zonal 
flow is being forced to a value well above that at which one of the waves 
excited by the topography is resonant. 

For a two-dimensional non-divergent flow between rigid boundaries the 
zonal momentum balance is 


iliil = [ v * ? *] + JL [ v *h*] - K ( [u] - [u] ) 
3t H 0 e 


( 1 ) 


Here brackets indicate zonal means and asterisks deviations from the zonal 
mean. The height of the topography is h and the mean depth of the fluid H 0 . 

The zonal flow is forced to a value [u] e in time k~ • f is the coriolis 

parameter and £ the relative vorticity. In CD’s argument a relative maximum of 

the second term on the r.h.s as a function of [u] , which is associated with a 

strong wave response near resonance, results in the second and third terms 
balancing each other in the ways just described. The first term on the r.h.s. 
plays no role. Although the argument is based on the linear wave response for 
a given zonal flow, CD goes on to show through numerical solutions that the full 
equations on a beta-plane channel also possess multiple stable states. 

The purpose of this paper is to present multiple states of the barotropic 
vorticity equation in which the balance is between the first and third terms 
on the r.h.s. of (1). Solutions of this type were also considered by CD in 
their discussion of thermally, rather than topographically, forces waves. 

Whereas in the case of topographic forcing the multiplicity arises from what 
they called "form-drag instability" in the wave-zonal flow interaction, in the 
"thermal forcing" case the associated instability appears to be the Rossby-wave 
instability discussed by Lorenz (1972), and the multiple states of the highly 
truncated model proved to be unstable when more degrees of freedom were added. 

The multiple statistically steady solutions described below are thus 
novel in that they do not involve form-drag instability, they have stable 
statistics in calculations with a large number of degrees of freedoms, and they 
occur on the sphere, with no artificial confinement in a resonant cavity. 

Further the solutions are obtained with no external forcing of the zonal flow. 

We will use the full non-linear equations for two-dimensional non-divergent 
motion between smooth, rigid boundaries on a sphere. The zonally asymmetric 
part of the flow is governed by 
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— I — C_i (uq - [uq]) + ~JL cos (J> (vq - [vq])} + F*, (2) 

acos(() 9X 9<{) 


* 

3C - - 
3t 


where q = C + f is the absolute vorticlty. We will consider only 
simple forcings of the form 


F* = k (Cg- ?*) , 


(3) 


where k is constant arid ? e a prescribed function of latitude and longitude, 
with 


[Cel =0 • 

With no topography (h* = 0), and no external forcing of the zonal flow ( [u ] e 
= 0) (1) becomes 


iinl = [v*s*] - k [u] , (4) 

3t 

a simple balance between frictional and eddy stresses, Eqs. (1) and (4) are 
solved numerically using finite differences on latitude-longitude coordinates. 

A resolution of 2° latitude by 2.5° longitude is used. Integrations are carried 
forward until a meaningful statistical steady state can be obtained, usually 
400 days. The forcing used is of the form 

2 

c = e exp - IAJlAo } sin(mX) . (6) 

e ° A<|> 

In all cases A<f> = 15°, <f> Q = 40°, and k = 1/15 d~^. 

Experiments with m = 3 . 

Experiments were performed with various wavenumbers and frequencies of forcing. 
Here we show results only for a stationary forcing in wavenumber 3, which illustrate 
the multiple states we found. Similar multiplicities occur for shorter waves with 
eastward propagating forcing. 

Figure 1 shows the zonal flows obtained for four values of £o when the calcu- 
lations were started from rest. The interpretation of these is straight forward: 
westerlies appear at the latitudes where the waves are forced and friction on 
the zonal wind balances eddy convergence of westerly momentum. Away from the 
forcing the forced wave is damped, depositing easterly momentum. This produces 
a zonal flow with critical latitudes for the forced stationary wave just north 
and south of the region of forcing. Note that easterlies are stronger on the 
poleward than on the equatorward side, as one would expect if the waves tend to 
transport equal angular momentum in both directions. 

As the forcing is increased this zonal flow pattern is amplified. However, 
since the zonal mean absolute vorticity gradient is being modified as the zonal 
flow changes, this amplification cannot continue indefinitely with increased 
forcing without a substantial change in the flow. Figure 2 shows the quantity 3m, 
as defined by Hoskins and Karoly (1981), 
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= — cos 2 .)) - £ 21 ' i> _i _i_ _1 (cos.)) [u] 
a a 9<|> cos <p 9<j> 

which is cos4> turns the absolute vorticity gradient. In the westerly jet, it is 
increased above its planetary value, and in the easterly jets it is depressed. 

For the largest forcing it is already slightly negative in the polar easterly 
jet, so that the inviscid stability condition is not satisfied. This is in fact 
the strongest forcing for which solutions like those in figure 1 were obtained. 
Above this value the solution is dramatically different, with no polar easterlies, 
a broad, strong westerly jet over the forcing, and balancing easterlies in the 
subtropics of forced hemisphere and mid-latitudes of the unforced hemisphere. 

This very sharp transition in the response suggested a bifurcation of the solution 
for which two stable branches may exist for some range of parameters. To explore 
this possibility the four experiments were repeated using initial conditions 
with a westerly jet of ~ 80 m/s at 40° in both hemispheres and easterlies in the 
tropics. The first three cases quickly evolved to the solutions shown in figure 
1. For the strongest forcing, however, a stable upper branch solution was 
obtained. Figure 3 shows the time evolution for this case when started from 
rest and when started from the strong westerly jet. The characteristics of the 
upper branch steady state just described are evident. 
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A New CO 2 Transmittance Parametization and its Impact on the GLA GCM 
R. L. Wobus, M.-L. C* Wu, and J. Susskind 


As reported last year (Wobus et al » , 1984), we are improving the Wu-Kaplan 
radiation parameterization (Krishnamurthy , 1982) used in the GLA GCM by replacing 
its fixed tables of CO 2 transmittance in the 15 ym band with models developed 
by regression on line-by-line transmittances. 

The transmittances between layers are modeled as products 

of effective sublayer transmittances (x£ k+i) : 


~ j 1 ~ , 

Ti,j \=i Tk,k+1 


Of the models tested offline, we select for the first sublayer (k=i) the 4-term 
model bilinearly dependent on layer temperature and surface pressure and for 
the remaining sublayers (k£L) the 2-term model linearly dependent on layer 
temperature. 


We have integrated the GLA GCM for 20 days starting at 0Z, January 21, 

1979, using the transmittance model described above. In the control run the 
fixed table of 15 ym CO 2 transmittances is used. The effect of the change 
of initial cooling rate is illustrated by (Fig.l) a map of the difference 
of 50 mb temperature after 6 hours. The cooling is reduced over high 
topography, where the fixed table underestimates the transmittance, and is 
reduced slightly throughout the tropics and the north polar area where 
the stratosphere is relatively cold. Over elevated topography the surface 
cooling (not shown) increases, also as expected. The final 10— day cross-section 
of temperature is shown in Fig. 2. The stratospheric temperature increases 
over a degree in the arctic and smaller amounts over Antarctica and elsewhere. 
Tropospheric equilibrium temperature response is obscured by time dependent 
differences in synoptic disturbances. 

The introduction of more accurate CO 2 transmittances allows us to evaluate 
two effects which are neglected in the present model. The increase of transmit- 
tance between levels as the layers grow thinner over elevated terrain affects 
both surface and stratospheric cooling rates. The change of CO 2 transmittance 
in the stratosphere with temperature, due to the inactivity of the hot bands at 
low temperatures, leads to warming over the winter pole. If this result is 
confirmed when the test runs are extended we plan to implement the new parameteri- 
zation in the GLA GCM. 
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50 mb temperature difference (new method-control) after 6 hours 
(6Z 21 January 1979), negative areas shaded, contours correspond 
to .l°K/day. 
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Fig. 2 50 mb temperature difference (new method-control), 10-day average 
after 10 days (6Z 31 January 1979 - 0Z 10 February 1979), negative 
areas shaded. 
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C. CLIMATE SENSITIVITY EXPERIMENTS 



A GCM STUDY OF THE ATMOSPHERIC RESPONSE 
TO TROPICAL SST ANOMALIES 

Max* J* Suarez 


The response of the atmosphere to anomalies in equatorial sea surface 
emperatures (SSTs) has received a great deal of attention in recent years, 
oth statistical studies using observational data and general circulation 
odel experiments have clearly established that the ocean temperature perturba- 
ions occurring during warm El Nino episodes cause a discernible climate 
ignal over much of the globe. A body of theory has also been developed that 
ttempts to explain the mid-latitude portion of the response in terms of global 
aleconnection patterns and propagating Rossby waves forced by the anomalous 
tmospheric heating. The most recent (1982/83) warm episode is of particular 
nterest because of the exceptional strength and highly unusual sequence of 
he SST anomalies. 

The purpose of the work presented here is to understand the evolution of 
he atmospheric anomalies associated with the most recent warm episode by the 
se of simulation studies with the UCLA general circulation model (GCM). Our 
pproach is to integrate the model using the observed sequence of SST anomalies 
uring 1982/1983, and compare it with a control run in which all boundary 
onditions vary from month to month as in the climatology. 

The control for the experiment was a three year simulation using seasonally 
arying climatological SSTs. The anomaly calculation was initialized from 15 
une of the first year of the control. From 15 June to 1 July, the SST was 
radually modified by the anomalies observed during June - July of 1982. 
rom 1 July on, the run was continued using the control’s SST plus the 1982-83 
nomalies. SSTs were varied daily, interpolating between monthly means. This 
econd integration was carried to the end of February of the second year ('83 
n the anomaly). Only anomalies over the tropical Pacific were used. 

Figures 1 and 2 summarize the model’s response. Seasonally averaged 
ST anomalies are shown together with the resulting precipitation and surface 
ind anomalies. During September - October - November (SON) maximum precipita- 
ion anomalies occur over the equator between the dateline and 120 °W with the 
reatest anomalies just east of the dateline. The precipitation anomaly is 
hus well west of the SST anomaly. By December-January-February (DJF) the 
esponse has shifted east some 30° and anomalies of ^ 6mm/day extend to the 
outh American coast. The total rainfall averaged over the tropical Pacific 
hanges little. Positive anomalies over the equator are compensated by nega- 
ive anomalies over a large horseshoe-shaped region poleward and westward of 
he increased rainfall. The major effect of the SST anomaly is thus a redis- 
ribution of the precipitation over the tropical Pacific. The anomaly in the 
urface wind, which must be doing the bulk of the anomalous transport of water 
apor, is also shown on the figure. Over the western half of the region of 
axiraum forcing - (180W-150W during SON and 150W-120W during DJF), - the result 
s qualitatively as one would expect from the simplest linear theory. Cyclonic 
urface flow north and south of the equator just west of the maximum forcing, 
nd strong westerly anomalies over the equator. East of the maximum rainfall 
nomaly, however, the response is very different from what one would expect 
rom linear theory. The model produces a strong meridional convergence with 
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particularly strong northerly flow just north of the equator. The linear 
response about a state of rest would lead to either flow away from the equator 
or, far east of the forcing, no meridional component. 



Surface Temperature Anomaly 


Figure 2 


Precipitation and 
Surface Wind Anomalies 







INFLUENCE OF LAND SURFACE PROCESSES ON THE INDIAN MONSOON 
- A NUMERICAL STUDY - 


Y. C. Sud and W. E. Smith* 

Laboratory for Atmospheres, Goddard Space Flight Center, Greenbelt, Md 20771 


Twelve July integrations are made with the 
GLAS GCM (Sud and Smith, 1984) to investigate the 
influence of local changes In the land surface 
fluxes, that may be produced by changes In the 
land surface vegetation, on the monsoon circula- 
tion and rainfall over the Indian subcontinent. 
The study consists of an ensemble mean of three 
integrations each for four separate cases: (a) 

Control Integrations with norraaly prescribed 
boundary conditions, called C; (b) three other 
integrations with the surface albedo increased 
(from 0.14 to 0.20) on the Indian subcontinent, 
called Ej; (c) another set of three integrations 
with the surface albedo Increased as above to- 
gether with a decrease of surface roughness 
(from 45 cm to 0.02 cm), called E 2 ; and (d) a 
last set of three integrations with higher sur- 
face albedo and low surface roughness as above, 
but with no evapotranspiration, called E 3 . Except 
for these changes in the anomaly region (delimit- 
ed by solid line, Fig. la) all other boundary 
conditions were prescribed, (Kalnay et al . , 1983); 
consequently, their time evaluation Is unaffected 
by the earth-atmosphere Interactions. 

Table 1 gives the regional summary of Impor- 
tant physical and diagnostic fields as well as 
the intrinsic model variability. Except for 
500 mb heights, all other differences are statis- 
tically significant. The most important fields 
are those of heat and moisture balance. Large 
differences were expected in these fields for E^ 
and E 3 ; nevertheless, the surface roughness ano- 
maly experiment E 2 has also shown significant 
differences in the surface radiation balance and 
the PBL fluxes. The rainfall differences for 
the three experiments are shown In Figs, la, b 
and c respectively. The surface albedo and sur- 
face roughness anomalies have produced comparable 
changes. For experiment E 2 , the most signifi- 
cant reductions in rainfall are found over north- 
northwest India (Fig. 2). The overall influence 
of the increased surface albedo is to reduce the 
rainfall over India; the rainfall differences E^ 
minus C show an average decrease of about 1.6 
mm/day with a maximum decrease of about 3 mm/day. 
This occurs due to relative subsidence (reduced 
rising motion) and diabatic cooling precisely in 
conformity with Charney's (1975) albedo hypothe- 
sis. For experiment E 2 the reduction in rainfall 
was even larger (2.7 mm/day) with a maximum 
decrease of about 5mm/day; it was caused by east- 
ward shift of rainfall maxima. 


*M/A-Com Sigma Data under contract No. NAS5-28074 




Fig. 1. Rainfall anomalies (Experiment-Control) 
mm/day for E^, E 2 and E 3 designated by 
a, b and c, negative contours are dashed. 
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The current findings could be easily under- 
stood by examining the moisture convergence 
fields for experiments Ej and Eo. In nature, the 
July rainfall anomalies over India are caused by 
changes in moisture convergence in the lower trop- 
osphere Often, small changes in the moisture 
flux (Vq) can produce a large change in the mois- 
ture convergence (V • Vq). This is indeed true 
for our simulations. We found that for smoother 
land the PBL wind vectors turned eastward which 
consequently produced the rainfall anomalies. 

TABLE 1. EXPT. MINUS CONTROL AND STD. DEVIATION 


Fields 


Expt. - Control 


St. Dev. 



E x - C 

E 2 - C 

E 3 - C 

|J L. . L/C V < 

Surf. Short Wave 
(ASW; W m -2 ) 

-17.9 

-18.9 

-12.3 

4.10 

Surf. Net Radn. 
(ANR; W m“ 2 ) 

-19.5 

-35.8 

-65.9 

7.92 

Surf. Sen. Flux 
(ASH; W m -2 ) 

-12.7 

-44.5 

45.5 

2.72 

Evapotranspiration 
(AE; mm/day) 

-0.24 

0.34 

-3.85 

0.10 

Percepitat ion 
(AP; mm/day) 

-1.62 

-2.7 

- 2.8 

0.68 

Surface Stress 
(At; N m” 2 ) 

-0.04 

-0.18 

-0.16 

0.03 

Vert Vel . 500 mb 
(Am; mb s“*) 

2.2 

3.0 

0.15 

1.14 

Diabatic heating 
(AQ; W m" 2 ) 

-57.2 

- 111.6 

12.1 

19.50 

Surface Temp. 
(AT S ; °C) 

-0.5 

- 0.8 

2.8 

0.17 

Ground Temp. 
(AT G ; °C) 

- 0.6 

2.1 

10.5 

0.11 

CeopoCentlal 500 mb 

-4.72 

-5.56 

5.72 

4.26 


(A<J>; GMP) 



The evapotranspiration is set to zero in 
experiment E 3 ; however, there is little change 
in the rainfall between experiments E 2 and E 3 
because the enhanced moisture convergence in the 
PBL, produced as a consequence of the increased 
sensible heating over dry land, compensated for 
the lack of evapotranspiration. This is an ex- 
ample of the typical 'heat island' effect, e.g. 
Malkus (1963). In July, the moisture for preci- 
pitation is advected from the nearby Indian 
Ocean; the increased sensible heating from non- 
evaporating land enhances this process. This 
is consistant with previous results (Shukla and 
Mintz, 1982; and Sud and Fennessy, 1984). 

These studies demonstrate that any modifica- 
tion of the biosphere in the Indian subcontinent 
may be expected to influence the monsoon circula- 
tion by altering the surface energy balance, ori- 
entation of motion fields and moisture conver- 
gence, and the hydrological cycle. The results 
further suggest that excessive land use via de- 
struction of vegetation would weaken the monsoon; 
conversely, increased surface roughness and more 
vegetation would strengthen the monsoon. One 
can also infer that even substantial changes in 
the biosphere may not Install a permanent desert 
over India, but a 3-5 mm/day reduction in the 
monsoon rainfall would be enough to seriously 
impact the agriculture. 

Despite various limitations of numerical 
models, we have confidence in the present results 
because they can be explained by a sequence of 
well understood physical mechanisms independent 
of the behavior of the GCM. We hypothesize that 
feedback effect of land surface processes may be 
partial contributors to the well documented pro- 
gression of the Thar desert into neighboring 
arable regions. The non— critical role of evapo— 
transporation for the monsoon rainfall Is another 
important result. 
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Fig. 2 


Rainfall anoraaly/model 's Variability 
plots: negative contours are dashed; 
significant regions are hatched. 
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INFLUENCE OF LAND-SURFACE ROUGHNESS ON ATMOSPHERIC CIRCULATION AND RAINFALL: 

A SENSITIVITY STUDY WITH A GCM 


Y. C. Sud, J. Shukla* and Y. Mintz* 

Laboratory for Atmospheres, Goddard Space Flight Center, Greenbelt, MD 20771 


1. INTRODUCTION - 

The dependence of atmospheric circulation on 
various aspects of land surface boundary condi- 
tions has been discussed for many years. However, 
it is only recently that numerical model simula- 
tion studies have been made to investigate Influ- 
ences of surface albedo and soil moisture on the 
circulation. Mintz (1984) has reviewed eleven 
recent studies on the subject. Subsequently, a 
dozen more studies have appeared in the litera- 
ture. All of them, without exception, have 
shown a strong influence of these parameters on 
the circulation and rainfall. In this study, we 
have chosen to investigate the sensitivity of 
atmospheric circulation to surface roughness of 
land. This work follows the recent study of Sud 
and Smith (1984b) on the Influence of surface- 
roughness of deserts. 

2. EXPERIMENTAL DESIGN 

An ensemble of three pairs (Control and Ex- 
periment) of runs were made with a GCM (general 
circulation model). The GCM employed for this 
study has been described by Randall (1982). It 
has been used for several studies in the past. 



Negative contours are dashed. 


♦Professor Department of Meteorology University 
of Maryland, College Park, Maryland 20742 



Fig. 2 Ensemble mean July moisture convergence 
differences: Experiment minus Control in 
mm d“^. Negative contours are dashed. 


However, the version of the GCM used in the pre- 
sent investigation Includes an ensemble PBL for- 
mulation by Sud and Smith (1984a) which In turn 
was derived from the Deardorff (1972) PBL para- 
meterization.. 

In the Control runs, the surface roughness 
values were: 45 cm for land, 0.02 cm for Ocean, 
and 0.01 cm for sea-ice. In the Experiment 
runs, the surface roughness over land was reduced 
to 0.02 cm which is a typical value for sandy 
deserts. For 1 km deep PBL this reduction in 
land surface roughness produces approximately a 
four fold reduction in the bulk aerodynamic drag 
coefficient. Such a range is typical of varia- 
tions in the biosphere. Apart from this one 
change all other Initial and boundary conditions 
were the same in each pair of runs. 

The two cases consist of integrations thru 
the month of July starting from two initial states 
from the NMC analysis of observations for June 15, 
1979, and June 15, 1980. The third initial state 
was for June 20, 1980; it was derived by averag- 
ing NMC analysis of observations from June 7 thru 
July 6, 1980. The differences in the July ensem- 
ble mean of the Controls and Experiments were 
analyzed. 
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L RESULTS 


Although there wa h little difference in the 
sea-Level pressure and 500 mb geopotential height, 
large differences were found in the low level 
winds, moisture convergence and rainfall. Large 
rainfall reductions are seen over the eastern 
U.S., Central America, and Amazons in South 
America, Fig. 1. The rainfall Increases over 
the CTCZ region in the raonsoonal North Africa 
and South India and China. Nevertheless, this 
increase is balanced by a corresponding reduction 
in rainfall to the north and south of the ITCZ. 
Whereas there was virtually no change in evapo- 
transpi_ ration, but the PBL moisture convergence 
(? • Vq) differences, Fig. 2, show some very 
remarkable correspondence with rainfall differ- 
ences. This one-to-one correspondence in the 
PBL moisture convergence and rainfall is not 
surprising because convective rain in July ob- 
tains most of the moisture from the PBL. This 
result suggests that it is the variation in deep 
cumulus convection which essentially responds to 
moisture convergence in the PBL. 

The vector plot of moisture flux (Vq) ident- 
ify the physical mechanism responsible for the 
change. It is found that reduced surface drag 
over smoother land suppresses the croas-isobaric 
moisture convergence into the low pressure regions 
at high latitudes, Fig. 3. On the other hand, low 
surface roughness increases the intensity of 
winds and moisture convergence into the ITCZ. 
This produces sinking and moisture divergence at 
the neighboring latitudes. Some large difference 
over south China Sea and western Pacific merely 
reflect a slight shift in the strong monsoonal 
circulation . 

The statistical significance of our results, 
represented by rainfall differences normalized 
with the natural variability of the model, is 
shown in Fig. 4. All the Americas are signifi- 
cantly affected. The rainfall in the monsoonal 
regions is also significantly affected. Over 
Europe and .Asia too, there are several regions 
of significant rainfall changes. 



Fig. 3 Ensemble mean July moisture flux vector 
differences: Experiment minus Control in 
(mm d~l)ra. Arrow lengths are discretised 
Small arrows: 0.5 x 10 6 to 5.0 x 10 6 
Medium arrows: 5.0 x 10 6 to 10.0 x 10 6 
Long arrows: greater than 10.0 x 10 6 



Fig. 4 Statistical significance of rainfall 
differences: (V precip/Vo model). 
Negative contours are dashed. 


Our study suggests that land-surface rough- 
ness is an important parameter. In nature, land 
surface roughness, surface albedo and evapotrans- 
piration are strongly influenced by vegetation. 
Therefore, it is vital that future simulations 
of weather and climate are made with GCMs that 
include a suitable parameterization of the vege- 
tation on land. Recently, Anthes (1984) has 
examined the possibility of enhancement of convec- 
tive precipitation by employing mesoscale varia- 
tion in vegetation. Baaed on our studies on the 
large scale, we infer that vegetation induced 
surface roughness too can support such a design. 
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III. SUMMER LECTURE SERIES 




The Structure and Dynamics of an Observed Moist Front 

by 

I. Orlanski 

Geophysical Fluid Dynamics Laboratory/NOAA 
Princeton University 
Princeton, New Jersey 08542 

The structure and dynamics of the moist cold front of 25-26 April 1979, 
the third observing day of the SESAME Experiment, are investigated through the 
use of a three-dimensional mesoscale numerical model. This work is one of the 
first studies in which model results are compared, in a one-to-one manner, 
with a detailed observational analysis, namely that of Ogura and Portis (1982) 
as taken from the SESAME observations. In addition, frontogenetical effects, 
both adiabatic and diabatic, are studied on a vertical cross-section through 
the front; similarities and differences with the adiabatic analysis of Ogura 
and Portis are discussed. 

The mesoscale model, which is described by Orlanski and Polinsky (1984), 
is initiated at 0000 GMT 25 April from the GFDL/FGGE global analysis and is 
integrated for 26 hours, using the same FGGE data at the synoptic times for 
boundary data, to produce the solution of 0200 GMT 26 April to be compared 
with the Ogura-Portis analysis. 

Many similarities exist between the modeled and analyzed fields, although 
the analysis tends to have weaker horizontal gradients due to the coarseness 
of the observational network. Vorticity and convergence near the surface are 
found to have the same magnitude in both the model solution and the analysis, 
in contrast to idealized frontogenesis models which predict vorticity to be 
much larger than convergence. The model produces a strong low-level jet ahead 
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of the front, in agreement with observations, when moist convection is active 
but with only a simple planetary boundary later used. Finally, vertical 
penetration of the frontal zone, indicated by the line of maximum vorticity in 
d vertical plane through the front, is close to the slope analyzed by Ogura 
and Portis. In fact, an analogous solution run without moisture produced much 
shallower penetration, suggesting that fluid parcels in the moist environment 
follow the more steeply inclined lines of constant equivalent potential tem- 
perature just as parcels in the dry case follow constant-e lines. 

A study was also made of the frontogenetical terms which act to enhance 
and weaken the horizontal potential temperature gradients within the front. 

The adiabatic terms, namely, the convergence, deformation, and tilting terms, 
show good agreement, with regard to sign and vertical structure, with those 
analyzed by Ogura and Portis. The low-level structure of each of the moist 
fields is close to that shown in the dry solution, with convergence and defor- 
mation being frontogenetical and tilting frontolytical near the surface. 
However, moist terms exhibit a much deeper penetration with the tilting field 
dominant in both the simulation and the analysis of observations. On the 
other hand, when diabatic and diffusive terms are obtained from the numerical 
solution, one finds that the tilting field is largely canceled by diabatic 
heating effects (Figure la). Hence, the complete frontogenetical function 
(Figure lb) consists of: frontogenesis ahead of the front at low levels, due 

to convergence and deformation, and at middle levels due to diabatic heating; 
and frontolysis at middle levels to the rear of the maximum vorticity line, 
due to tilting and diffusive effects. 
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(b) FRONTOGENETICAL FUNCTION 



Figure 1 — Vertical cross-sections through the simulated cold front of 26 April 1979 showing (a) 
combined diabatic and tilting terms and (b) total frontogenetical function, both in units of 10 -10 
K m _ *s - Frontogenetical effects (stippling) are shown at ground level and upper levels ahead of the 
front, whereas frontolytical effects (cross-hatching) are shown behind the front. The effects at the 
ground produced by adiabatic processes (convergence and deformation) and the frontolytical effects 
(tilting) had been speculated from idealized studies. The upper frontogenetical effects ahead of the 
front produced by latent heat release (diabatic processes) are shown for the first time by this study. 
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Some Recent Developments in Numerical Modelling at ECMWF 


A.J. Simmons 

European Centre for Medium Range Weather Forecasts, 
Reading , Berks . f U . K . 


Extended Abstract 


1 . Introduction 

A new atmospheric model was introduced into operational forecasting at 
ECMWF on 21 April 1983. The principal differences between this model and 
the Centre's first operational model were in the adiabatic formulation, 
which, in the new model ,' includes use of a spectral representation in the 
horizontal, a more general vertical coordinate, and a modified, more-efficient, 
time-stepping scheme. In addition, new programming techniques and standards 
were adopted to facilitate both the model's use as a research tool and 
its adaptation to make full use of the features of the recently-acquired 
CRAY X-MP computer. A number of revisions were also made in detailed aspects 
of the formulation of the parameterisation schemes. The operational change 
to this new model was accompanied by a second important change, namely 
the use of a higher 'envelope* orography in the lower boundary conditions 
of the model. An account of the design and early performance of the new 
model has been given by Simmons and Jarraud (1984). 

Since its operational implementation, effort has continued to be directed 
towards the development of the new model. This effort can be divided into 
two categories. The first concerns development directly related to problems 
in the operational performance of the model. The second concerns the develop- 
ment preparatory to the future operational implementation of a high resolution 
version of the model on the CRAY X-MP. 

2 . Development of the Numerical Formulation of the Operational Model 
(a) The representation of orography 

In the version of the new model originally introduced into operational 
forecasting, the spectrally-fitted envelope orography was modified iteratively 
to reduce the height of the model orography in coastal areas, and thereby 
increase the use of synoptic observations from such areas in the data assimilation 
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oubsequent tests showed that any beneficial impact of such data was more 
than masked by the impact of an unrealistic model orography in near-coastal 
regions. It was thus decided to use an orography defined by a direct spectral 
fit of the envelope orography. In the construction of the latter, the 
standard-deviation of the sub grid-scale orography was not used to enhance 
the mean at model "sea" points (points for which more than 50% of the 
surrounding grid-square is in reality covered by sea). Improvements in 
local synoptic features were evident near regions of significant change 
in the orography, and some of the precipitation biases were reduced. 

In addition to making the above operational changes, the overall performance 
of the envelope orography in routine prediction has been monitored. Comparing 
the results of objective verifications for the summers of 1982 and 1983 
showed the improvement in quality in surface forecasts anticipated from 
previous comparisons of grid-point and spectral forecasts. However, a 
clear improvement at 500 mb of the type found previously to result from 
use of the envelope orography in winter forecasts was not evident. Although 
results for 1984 are substantially better and suggest that such interannual 
comparisons can be misleading due to variations in predictability in different 
synoptic situations, experimentation has indicated a small detrimental 
effect of the envelope orography in some summer forecasts. Some further 
investigation has been carried out, but it was decided to concentrate 
attention on co-ordinating this work with the experimental programme for 
the development of a higher resolution operational model. 

(b) The semi-implicit gravity-wave treatment and th 5 _h°£ iz 25 tal 
diffusion of divergence 

Noise generated in the model vertical velocity field in conjunction with 
intense convection has been found to be a recurrent feature of forecasts, 
particularly in the short-range. In addition, noise has occasionally been 
found to develop in the data assimilation in association with anomalously 
strong mid-latitude jet streams, and a case of computational instability 
was encountered operationally when the problem was exacerbated by serious 
errors in satellite temperature data. Tests revealed sensitivity to the 
treatment of gravity waves by the semi-implicit time scheme. Pending a 
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fuller investigation of this problem, a solution has been adopted which 
involves a minor modification of the semi-implicit scheme, together with 
an enhanced (by a factor of 10) coefficient of horizontal diffusion for 
the divergence field. The position of this diffusion within the model 
code was also corrected so that it was executed after completion of the 
semi-implicit correction of the temperature and surface pressure. The 
overall impact on the accuracy of forecasts was extremely small, although 
some indication of excessive damping of tropical features has been found. 

c) Other operational changes 

Other minor changes in numerical formulation have been introduced in the 
light of operational experience. Horizontal diffusion is now applied not 
to temperature itself, but to the deviation of temperature from a reference 
profile in such a way as to approximate diffusion on surfaces of constant 
pressure. This acts to suppress excessive precipitation over mountainous 
areas in summer forecasts, which was formerly triggered by an erroneous 
warming of mountain tops associated with the diffusion of temperature along 
the terrain-following coordinate surfaces of the model. A small increase 
in the time-filter coefficient and a modification of the reference wind profile 
used for the semi-implicit treatment of vorticity and moisture advection 
have enhanced the computational stability of the model. These changes 
have had negligible impact on forecast quality. 

3. Experimentation at High Resolution 

Operational forecasting is currently carried out using triangular truncation 
at wavenumber 63 (T63) in the horizontal and a 16-level vertical representation 
The acquisition of the CRAY X-MP computer allows use of a resolution up 
to about T106 using 16 levels, or about T95 if vertical resolution were 
to be increased to 20 levels. An experimental programme to assess the 
performance of the forecast model at higher horizontal and vertical resolution, 
and the sensitivity to orographic representation at such resolutions, has 
thus begun. 

First results indicate little sensitivity to an increase in vertical resolution 
but an overall improvement from an increase in horizontal resolution, at 
least in the extratropical Northern Hemisphere. Differences in conventional 
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objective verifications are not large in general but as no attempt was 
made m these first experiments to optimize initial and boundary conditions 
or model parameters, more substantial gains can be expected for the future. 
Improvement is more evident for the medium and smaller scales of motion, 
and can be clearly seen in the few precipitation forecasts that have been 
examined to date. Use of an envelope orography appears more beneficial 
at higher horizontal resolution, the small detrimental impact found in 
summer T63 forecasts disappearing at T106 resolution in the limited number 
of cases that have been completed. Further experimentation planned for 
the coming months should help quantify these conclusions and aid understanding 
of the sensitivity to orography. 
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LARGE-EDDY SIMULATION 
IN 

PLANETARY BOUNDARY-LAYER RESEARCH 
3. C. Wyngaard 

National Center for Atmospheric Research 
Boulder, CO 80307 

We do planetary boundary-layer research for the same reason we do 
research on any turbulent flow: to develop models (i.e., parameteriza- 

tions) of its structure and dynamics. While the traditional approach is 
through experiments, two developments over the past decade bear on the 
continuance of this tradition. 

The first is that progress in PBL research seems to have reached a 
plateau, in the view of a WMO working group (Andrd et al., 1982). Ex- 
periments are very expensive (the 1968 Kansas program required three 
previous trials, 15 persons in the field, about 10 persons over the next 
2-3 years in data analysis — all for measurements in the first 32 m!) and 
increasingly difficult to justify. Furthermore, inherent uncertainty-- 
the inevitable difference between time-averaged and ensemble-averaged 
structure — becomes much more important higher in the PBL, causing severe 
data scatter problems. As a result, when we look at new experimental 
challenges, such as reliably parameterizing subgrid-scale transport 
through a cloud layer, we can easily become discouraged; such problems 
seem out of experimental reach! 

The second development is the growing realization that large-eddy 
simulation, or LES, can provide data bases on PBL structure, and thereby 
supplement and extend the traditional direct measurements. LES refers 
to three-dimensional, fine-mesh, time-dependent turbulent flow simula- 
tion, with the spatial grid fine enough to resolve the energy-containing 
eddies. Deardorff (1970) demonstrated its potential first in turbulent 
channel flow, and then began his series of PBL simulations which pro- 
duced, for example, a study of the evolution of day 33 in the Wangara 
experiment (Deardorff, 1974). LES is used rather widely today in fluid 
mechanics research problems ranging from engineering flows (Rogallo and 
Moin, 1984) to severe storms (Klemp and Rotunno, 1983). A recent 
working-group report (Wyngaard, 1984) examines LES in detail and recom- 
mends its wider use as a tool in PBL research. 

In the Mesoscale Research Section of NCAR we are currently using 
LES to study the structure and dynamics of the convective boundary 
layer. In one study (Wyngaard and Brost, 1984) we simulated the 
vertical transport of a conservative, passive scalar. Because of the 
(mathematical) linearity of the governing transport equation for the 
scalar, we were able to represent the scalar field as the sum of "top- 
down" and "bottom-up" components. We isolated the statistical proper- 
ties of these two processes and found them to be significantly differ- 
ent; "bottom-up" diffusion has a larger (by a factor of 2 - 4) eddy 
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diffusivity. We attributed this difference to the asymmetry in the 
vertical profile of buoyant production of turbulence, and to the differ- 
ence in boundary conditions at bottom and top. 

In a second study (Moeng and Wyngaard, 1984) we studied in more 
detail the statistics of top-down and bottom-up scalar fields. We found 
substantial differences between them, due again, presumably, to the 
asymmetry in the convective PBL. We developed a generalization of 
mixed-layer scaling which allows one to include the effects of top-down 
diffusion (i.e., the effects of entrainment). 

Chin-Hoh Moeng is now extending this work to the study of cloud- 
topped boundary layers. We believe that LES has high potential for 
application to a very wide range of PBL problems, including turbulent 
dispersion, the stable PBL, and studies of mean wind profiles in baro- 
clinic cases. Some of tomorrow's field programs might take place on the 
supercomputer! 
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